Reason for Selection Protected natural areas in urban environments provide urban residents a nearby place to connect with nature and offer refugia for some species. They help foster a conservation ethic by providing opportunities for people to connect with nature, and also support ecosystem services like offsetting heat island effects (Greene and Millward 2017, Simpson 1998), water filtration, stormwater retention, and more (Hoover and Hopton 2019). In addition, parks, greenspace, and greenways can help improve physical and psychological health in communities (Gies 2006). Urban park size complements the equitable access to potential parks indicator by capturing the value of existing parks.Input DataSoutheast Blueprint 2024 extentFWS National Realty Tracts, accessed 12-13-2023Protected Areas Database of the United States(PAD-US):PAD-US 3.0 national geodatabase -Combined Proclamation Marine Fee Designation Easement, accessed 12-6-20232020 Census Urban Areas from the Census Bureau’s urban-rural classification; download the data, read more about how urban areas were redefined following the 2020 censusOpenStreetMap data “multipolygons” layer, accessed 12-5-2023A polygon from this dataset is considered a beach if the value in the “natural” tag attribute is “beach”. Data for coastal states (VA, NC, SC, GA, FL, AL, MS, LA, TX) were downloaded in .pbf format and translated to an ESRI shapefile using R code. OpenStreetMap® is open data, licensed under theOpen Data Commons Open Database License (ODbL) by theOpenStreetMap Foundation (OSMF). Additional credit to OSM contributors. Read more onthe OSM copyright page.2021 National Land Cover Database (NLCD): Percentdevelopedimperviousness2023NOAA coastal relief model: volumes 2 (Southeast Atlantic), 3 (Florida and East Gulf of America), 4 (Central Gulf of America), and 5 (Western Gulf of America), accessed 3-27-2024Mapping StepsCreate a seamless vector layer to constrain the extent of the urban park size indicator to inland and nearshore marine areas <10 m in depth. The deep offshore areas of marine parks do not meet the intent of this indicator to capture nearby opportunities for urban residents to connect with nature. Shallow areas are more accessible for recreational activities like snorkeling, which typically has a maximum recommended depth of 12-15 meters. This step mirrors the approach taken in the Caribbean version of this indicator.Merge all coastal relief model rasters (.nc format) together using QGIS “create virtual raster”.Save merged raster to .tif and import into ArcPro.Reclassify the NOAA coastal relief model data to assign areas with an elevation of land to -10 m a value of 1. Assign all other areas (deep marine) a value of 0.Convert the raster produced above to vector using the “RasterToPolygon” tool.Clip to 2024 subregions using “Pairwise Clip” tool.Break apart multipart polygons using “Multipart to single parts” tool.Hand-edit to remove deep marine polygon.Dissolve the resulting data layer.This produces a seamless polygon defining land and shallow marine areas.Clip the Census urban area layer to the bounding box of NoData surrounding the extent of Southeast Blueprint 2024.Clip PAD-US 3.0 to the bounding box of NoData surrounding the extent of Southeast Blueprint 2024.Remove the following areas from PAD-US 3.0, which are outside the scope of this indicator to represent parks:All School Trust Lands in Oklahoma and Mississippi (Loc Des = “School Lands” or “School Trust Lands”). These extensive lands are leased out and are not open to the public.All tribal and military lands (“Des_Tp” = "TRIBL" or “Des_Tp” = "MIL"). Generally, these lands are not intended for public recreational use.All BOEM marine lease blocks (“Own_Name” = "BOEM"). These Outer Continental Shelf lease blocks do not represent actively protected marine parks, but serve as the “legal definition for BOEM offshore boundary coordinates...for leasing and administrative purposes” (BOEM).All lands designated as “proclamation” (“Des_Tp” = "PROC"). These typically represent the approved boundary of public lands, within which land protection is authorized to occur, but not all lands within the proclamation boundary are necessarily currently in a conserved status.Retain only selected attribute fields from PAD-US to get rid of irrelevant attributes.Merged the filtered PAD-US layer produced above with the OSM beaches and FWS National Realty Tracts to produce a combined protected areas dataset.The resulting merged data layer contains overlapping polygons. To remove overlapping polygons, use the Dissolve function.Clip the resulting data layer to the inland and nearshore extent.Process all multipart polygons (e.g., separate parcels within a National Wildlife Refuge) to single parts (referred to in Arc software as an “explode”).Select all polygons that intersect the Census urban extent within 0.5 miles. We chose 0.5 miles to represent a reasonable walking distance based on input and feedback from park access experts. Assuming a moderate intensity walking pace of 3 miles per hour, as defined by the U.S. Department of Health and Human Service’s physical activity guidelines, the 0.5 mi distance also corresponds to the 10-minute walk threshold used in the equitable access to potential parks indicator.Dissolve all the park polygons that were selected in the previous step.Process all multipart polygons to single parts (“explode”) again.Add a unique ID to the selected parks. This value will be used in a later step to join the parks to their buffers.Create a 0.5 mi (805 m) buffer ring around each park using the multiring plugin in QGIS. Ensure that “dissolve buffers” is disabled so that a single 0.5 mi buffer is created for each park.Assess the amount of overlap between the buffered park and the Census urban area using “overlap analysis”. This step is necessary to identify parks that do not intersect the urban area, but which lie within an urban matrix (e.g., Umstead Park in Raleigh, NC and Davidson-Arabia Mountain Nature Preserve in Atlanta, GA). This step creates a table that is joined back to the park polygons using the UniqueID.Remove parks that had ≤10% overlap with the urban areas when buffered. This excludes mostly non-urban parks that do not meet the intent of this indicator to capture parks that provide nearby access for urban residents. Note: The 10% threshold is a judgement call based on testing which known urban parks and urban National Wildlife Refuges are captured at different overlap cutoffs and is intended to be as inclusive as possible.Calculate the GIS acres of each remaining park unit using the Add Geometry Attributes function.Buffer the selected parks by 15 m. Buffering prevents very small and narrow parks from being left out of the indicator when the polygons are converted to raster.Reclassify the parks based on their area into the 7 classes seen in the final indicator values below. These thresholds were informed by park classification guidelines from the National Recreation and Park Association, which classify neighborhood parks as 5-10 acres, community parks as 30-50 acres, and large urban parks as optimally 75+ acres (Mertes and Hall 1995).Assess the impervious surface composition of each park using the NLCD 2021 impervious layer and the Zonal Statistics “MEAN” function. Retain only the mean percent impervious value for each park.Extract only parks with a mean impervious pixel value <80%. This step excludes parks that do not meet the intent of the indicator to capture opportunities to connect with nature and offer refugia for species (e.g., the Superdome in New Orleans, LA, the Astrodome in Houston, TX, and City Plaza in Raleigh, NC).Extract again to the inland and nearshore extent.Export the final vector file to a shapefile and import to ArcGIS Pro.Convert the resulting polygons to raster using the ArcPy Feature to Raster function and the area class field.Assign a value of 0 to all other pixels in the Southeast Blueprint 2024 extent not already identified as an urban park in the mapping steps above. Zero values are intended to help users better understand the extent of this indicator and make it perform better in online tools.Use the land and shallow marine layer and “extract by mask” tool to save the final version of this indicator.Add color and legend to raster attribute table.As a final step, clip to the spatial extent of Southeast Blueprint 2024.Note: For more details on the mapping steps, code used to create this layer is available in theSoutheast Blueprint Data Downloadunder > 6_Code. Final indicator valuesIndicator values are assigned as follows:6= 75+ acre urban park5= 50 to <75 acre urban park4= 30 to <50 acre urban park3= 10 to <30 acre urban park2=5 to <10acreurbanpark1 = <5 acre urban park0 = Not identified as an urban parkKnown IssuesThis indicator does not include park amenities that influence how well the park serves people and should not be the only tool used for parks and recreation planning. Park standards should be determined at a local level to account for various community issues, values, needs, and available resources.This indicator includes some protected areas that are not open to the public and not typically thought of as “parks”, like mitigation lands, private easements, and private golf courses. While we experimented with excluding them using the public access attribute in PAD, due to numerous inaccuracies, this inadvertently removed protected lands that are known to be publicly accessible. As a result, we erred on the side of including the non-publicly accessible lands.The NLCD percent impervious layer contains classification inaccuracies. As a result, this indicator may exclude parks that are mostly natural because they are misclassified as mostly impervious. Conversely, this indicator may include parks that are mostly impervious because they are misclassified as mostly

Survey participants plotted activity points using an interactive mapping tool.The 2012 Northeast Recreational Boater Survey was conducted by SeaPlan, the Northeast Regional Ocean Council (NROC), states’ coastal agencies, marine trade associations composed of many private industry representatives, and the First Coast Guard District. The methodology for the 2012 Northeast Recreational Boater Survey follows a protocol similar to the 2010 Massachusetts Survey with modifications based on the lessons learned and recommendations suggested in the Massachusetts Survey Final Report.The methodology consists of surveying a random sample of selected boat owners throughout the Northeast through a series of monthly online surveys. The surveying period lasted throughout the 2012 boating season (May 1 through October 31, 2012), which was identified by the advisory committee (consisting of NROC and representatives from the recreational boating industry).The project team decided to use a random sample survey approach because it successfully gathered statistically robust economic and spatial data on recreational boating activity by Massachusetts registered boaters during the 2010 boating season. This was also the only approach that would allow for the calculation of statistically robust economic impact estimates for both the states and the region, which was identified as a priority (along withspatial data) by both NROC and the boating industry.

This dataset can be used by coastal planners in ocean planning activities to develop a better understanding of how and where humans use the ocean in the Northeast to inform regional ocean planning and minimize ocean use conflicts. This effort also fulfilled a recommendation from the 2010 Massachusetts Survey to expand the survey’s geographic range to the Northeast Region, allowing for the capture of interstate traffic between states in the Northeast. Furthermore, this dataset can also be used by the boating industry to show the importance of recreational boating to the region and to inform business planning.

Supplemental Information; SURVEY SAMPLING METHODOLOGY - The sample for this survey came from seven databases, including the U.S. Coast Guard Documented Vessel Database and databases of state registered boaters from New York, Connecticut, Rhode Island, Massachusetts, New Hampshire, and Maine. Recreational boaters who owned vessels that met the following criteria were eligible for the survey: * Registration: Currently registered with a state in the Northeast and/or registered as a documented vessel with the U.S. Coast Guard, with a hailing port in the Northeast * Primary Use: Recreational use designation * Length: At least 10 feet in length * Saltwater (if specified; only Maine and New Hampshire required this information) * Location: Located in a “coastal county”. The survey team defined “coastal counties” as those that border saltwater, or those that were highlighted by state coastal planners as likely containing large amount of saltwater boating activity. Based on the 2010 Massachusetts Survey and budgetary considerations, the project team determined an overall sample size that would provide sufficient spatial and economic data for both each state, as well as the whole Northeast. Because of the, at times, large discrepancies between the number of eligible boats in some states, the team decided that certain states with fewer eligible boats should also have a supplemental sample of boats in addition to the pure random sample. To ensure the sample represented the total population of registered boats in the Northeast, the sampling method included considerations of state, geography and size class. Of the 373,766 boats eligible for the survey, the base of randomly sampled boats included 50,000 boats from across all six states. In addition to this base, the survey team sampled 17,772 boats as a supplemental sample, including: 1,772 boats of 26 feet in length or more from across all six states to increase the number of large boats in the sample, and 16,000 additional boats to ensure each state had enough responses for the statistical analysis. These included 10,000 boats from Maine, 2,500 boats from Rhode Island, 2,000 boats from New Hampshire and 1,500 boats from Connecticut. This resulted in a total of 67,772 boaters invited to participate in the study. Boater Recruitment and Response: In the survey invitation package, the survey team also sent invited boaters a questionnaire to verify eligibility to participate in the survey. Eligibility requirements consist of: boat is used in saltwater; boat is used for recreational purposes; and boaters have access to the internet with a working email address. 12,218 boaters responded to the invitation; however only 7,800 of these respondents were found to meet all of the above criteria. From this sample, 4,297 individual boaters completed at least one monthly survey. Surveying Process: The study consisted of six monthly surveys and one end of season survey. The online monthly surveys gathered spatial and economic data on recreational boating activity that occurred during the previous month. The online survey had two parts: 1) a survey with questions about general boating activity during the previous month, and the boater’s last trip of the month (specifically focusing on spending), and 2) a mapping application developed by Ecotrust where boaters plotted their boating route and identified any areas where they participated in activities, such as fishing, diving, wildlife viewing, swimming and relaxing at anchor. The end of season survey gathered a variety of information that could not be gathered in the monthly surveys. The end of season survey contained questions about yearly boating-related expenditures (e.g., dockage, storage, taxes, yearly maintenance), feedback on the survey itself, and general boating-related questions (e.g. whether boaters have taken a boating safety course). Density Analysis: The density analysis described in the following paragraphs was vetted by a technical advisory team consisting of representatives from the Massachusetts Office of Coastal Zone Management (MA CZM), NROC, Maine Coastal Program and Applied Science Associates (ASA) and was based on mapping and analysis protocols from the 2010 Massachusetts Survey. To develop the density layer, vessel routes were drawn in WGS 1984 in the Ecotrust mapping application and were imported into Excel, then ArcMap using a data frame in that coordinate system. Routes from the random sample were selected from that data layer, and the data layer was re-projected into two separate shapefiles, one in UTM 18 and one in UTM 19. A line density analysis using a 250 m square grid cell with a 675 m neighborhood was applied to each shapefile. The 675 m neighborhood was applied to account for inherent user error in the mapping tool. The line density analysis resulted in a raster grid for each UTM zone. Each raster was clipped by the boundaries of its UTM zone, re-projected into the North American Albers Equal Area Conic Projection, and the separate rasters were mosaicked together. At the boundary of the two raster grids there was a line of cells with no data value. This was a result of mosaicking rasters that originated in different coordinate systems. To approximate values in the blank cells, each blank cell was populated by a value from a focal statistics calculation. The focal statistics expression took the mean of all cells in a 4x4 neighborhood around each blank cell. The values were then converted to Z-scores using the raster calculator by taking the log of the density values, subtracting the mean value, and dividing the resulting value by the standard deviation of the value. This layer was clipped again using the NOAA medium resolution shoreline dataset. DATA PROCESSING Processing environment: ArcGIS 10.05, Windows 7 Ultimate SP5, Intel Xeon CPU Process Steps Description 1 Raw routes from mapping application imported into ArcMap 2 Routes from random sample selected using select by attributes query 3 Routes projected into two separate shapefiles (UTM Zones 18 & 19) 4 LINE DENSITY tool in spatial analyst applied to each shapefile using a 250 m square grid with a 675 m neighborhood 5 Resulting rasters clipped to their respective UTM Zones using the EXTRACT BY MASK tool 6 Rasters reprojected to North America Albers Equal Area Conic Projection, using PROJECT tool 7 MOSAIC tool used to merge rasters 8 Focal mean expression (4x4 neighborhood) used to approximate and fill cells with no data at the boundary between mosaicked rasters 9 Raster calculator used to calculated Z-scores ([(Ln(Value))-Mean]/Std. Deviation) 10 Raster clipped by NOAA Medium Resolution Shoreline data using EXTRACT BY POLYGON tool QUALITY PROCESS Attribute Accuracy: The lines used to generate the density grid were derived from a mapping tool used by boaters to reconstruct their boating routes. To ensure that boaters included their round-trip route the mapping applications would send the user an error message asking them to re-plot the route or the program would automatically return the route to the starting point. This application also restricted the scale at which users could draw their routes, reducing the amount of error that could occur from plotting routes at too small a scale. Clipping this layer with a regional ocean shapefile derived from the NOAA medium resolution shoreline dataset excluded route density resulting from routes drawn over land, in freshwater, or outside of northeastern waters. Logical Consistency: None Completeness: Only reported routes from the random sample were included. Routes from the supplemental sample were excluded from this analysis. Route density occurring over land, freshwater areas, or outside northeastern waters was excluded by the final geoprocessing step. Positional Accuracy: The positional accuracy of the routes is dependent on the individual reporting routes through the

Overview

This dataset is the repository for the following paper submitted to Data in Brief:

Kempf, M. A dataset to model Levantine landcover and land-use change connected to climate change, the Arab Spring and COVID-19. Data in Brief (submitted: December 2023).

The Data in Brief article contains the supplement information and is the related data paper to:

Kempf, M. Climate change, the Arab Spring, and COVID-19 - Impacts on landcover transformations in the Levant. Journal of Arid Environments (revision submitted: December 2023).

Description/abstract

The Levant region is highly vulnerable to climate change, experiencing prolonged heat waves that have led to societal crises and population displacement. Since 2010, the area has been marked by socio-political turmoil, including the Syrian civil war and currently the escalation of the so-called Israeli-Palestinian Conflict, which strained neighbouring countries like Jordan due to the influx of Syrian refugees and increases population vulnerability to governmental decision-making. Jordan, in particular, has seen rapid population growth and significant changes in land-use and infrastructure, leading to over-exploitation of the landscape through irrigation and construction. This dataset uses climate data, satellite imagery, and land cover information to illustrate the substantial increase in construction activity and highlights the intricate relationship between climate change predictions and current socio-political developments in the Levant.

Folder structure

The main folder after download contains all data, in which the following subfolders are stored are stored as zipped files:

“code” stores the above described 9 code chunks to read, extract, process, analyse, and visualize the data.

“MODIS_merged” contains the 16-days, 250 m resolution NDVI imagery merged from three tiles (h20v05, h21v05, h21v06) and cropped to the study area, n=510, covering January 2001 to December 2022 and including January and February 2023.

“mask” contains a single shapefile, which is the merged product of administrative boundaries, including Jordan, Lebanon, Israel, Syria, and Palestine (“MERGED_LEVANT.shp”).

“yield_productivity” contains .csv files of yield information for all countries listed above.

“population” contains two files with the same name but different format. The .csv file is for processing and plotting in R. The .ods file is for enhanced visualization of population dynamics in the Levant (Socio_cultural_political_development_database_FAO2023.ods).

“GLDAS” stores the raw data of the NASA Global Land Data Assimilation System datasets that can be read, extracted (variable name), and processed using code “8_GLDAS_read_extract_trend” from the respective folder. One folder contains data from 1975-2022 and a second the additional January and February 2023 data.

“built_up” contains the landcover and built-up change data from 1975 to 2022. This folder is subdivided into two subfolder which contain the raw data and the already processed data. “raw_data” contains the unprocessed datasets and “derived_data” stores the cropped built_up datasets at 5 year intervals, e.g., “Levant_built_up_1975.tif”.

Code structure

1_MODIS_NDVI_hdf_file_extraction.R

This is the first code chunk that refers to the extraction of MODIS data from .hdf file format. The following packages must be installed and the raw data must be downloaded using a simple mass downloader, e.g., from google chrome. Packages: terra. Download MODIS data from after registration from: https://lpdaac.usgs.gov/products/mod13q1v061/ or https://search.earthdata.nasa.gov/search (MODIS/Terra Vegetation Indices 16-Day L3 Global 250m SIN Grid V061, last accessed, 09th of October 2023). The code reads a list of files, extracts the NDVI, and saves each file to a single .tif-file with the indication “NDVI”. Because the study area is quite large, we have to load three different (spatially) time series and merge them later. Note that the time series are temporally consistent.

2_MERGE_MODIS_tiles.R

In this code, we load and merge the three different stacks to produce large and consistent time series of NDVI imagery across the study area. We further use the package gtools to load the files in (1, 2, 3, 4, 5, 6, etc.). Here, we have three stacks from which we merge the first two (stack 1, stack 2) and store them. We then merge this stack with stack 3. We produce single files named NDVI_final_*consecutivenumber*.tif. Before saving the final output of single merged files, create a folder called “merged” and set the working directory to this folder, e.g., setwd("your directory_MODIS/merged").

3_CROP_MODIS_merged_tiles.R

Now we want to crop the derived MODIS tiles to our study area. We are using a mask, which is provided as .shp file in the repository, named "MERGED_LEVANT.shp". We load the merged .tif files and crop the stack with the vector. Saving to individual files, we name them “NDVI_merged_clip_*consecutivenumber*.tif. We now produced single cropped NDVI time series data from MODIS.
The repository provides the already clipped and merged NDVI datasets.

4_TREND_analysis_NDVI.R

Now, we want to perform trend analysis from the derived data. The data we load is tricky as it contains 16-days return period across a year for the period of 22 years. Growing season sums contain MAM (March-May), JJA (June-August), and SON (September-November). December is represented as a single file, which means that the period DJF (December-February) is represented by 5 images instead of 6. For the last DJF period (December 2022), the data from January and February 2023 can be added. The code selects the respective images from the stack, depending on which period is under consideration. From these stacks, individual annually resolved growing season sums are generated and the slope is calculated. We can then extract the p-values of the trend and characterize all values with high confidence level (0.05). Using the ggplot2 package and the melt function from reshape2 package, we can create a plot of the reclassified NDVI trends together with a local smoother (LOESS) of value 0.3.
To increase comparability and understand the amplitude of the trends, z-scores were calculated and plotted, which show the deviation of the values from the mean. This has been done for the NDVI values as well as the GLDAS climate variables as a normalization technique.

5_BUILT_UP_change_raster.R

Let us look at the landcover changes now. We are working with the terra package and get raster data from here: https://ghsl.jrc.ec.europa.eu/download.php?ds=bu (last accessed 03. March 2023, 100 m resolution, global coverage). Here, one can download the temporal coverage that is aimed for and reclassify it using the code after cropping to the individual study area. Here, I summed up different raster to characterize the built-up change in continuous values between 1975 and 2022.

6_POPULATION_numbers_plot.R

For this plot, one needs to load the .csv-file “Socio_cultural_political_development_database_FAO2023.csv” from the repository. The ggplot script provided produces the desired plot with all countries under consideration.

7_YIELD_plot.R

In this section, we are using the country productivity from the supplement in the repository “yield_productivity” (e.g., "Jordan_yield.csv". Each of the single country yield datasets is plotted in a ggplot and combined using the patchwork package in R.

8_GLDAS_read_extract_trend

The last code provides the basis for the trend analysis of the climate variables used in the paper. The raw data can be accessed https://disc.gsfc.nasa.gov/datasets?keywords=GLDAS%20Noah%20Land%20Surface%20Model%20L4%20monthly&page=1 (last accessed 9th of October 2023). The raw data comes in .nc file format and various variables can be extracted using the [“^a variable name”] command from the spatraster collection. Each time you run the code, this variable name must be adjusted to meet the requirements for the variables (see this link for abbreviations: https://disc.gsfc.nasa.gov/datasets/GLDAS_CLSM025_D_2.0/summary, last accessed 09th of October 2023; or the respective code chunk when reading a .nc file with the ncdf4 package in R) or run print(nc) from the code or use names(the spatraster collection).
Choosing one variable, the code uses the MERGED_LEVANT.shp mask from the repository to crop and mask the data to the outline of the study area.
From the processed data, trend analysis are conducted and z-scores were calculated following the code described above. However, annual trends require the frequency of the time series analysis to be set to value = 12. Regarding, e.g., rainfall, which is measured as annual sums and not means, the chunk r.sum=r.sum/12 has to be removed or set to r.sum=r.sum/1 to avoid calculating annual mean values (see other variables). Seasonal subset can be calculated as described in the code. Here, 3-month subsets were chosen for growing seasons, e.g. March-May (MAM), June-July (JJA), September-November (SON), and DJF (December-February, including Jan/Feb of the consecutive year).
From the data, mean values of 48 consecutive years are calculated and trend analysis are performed as describe above. In the same way, p-values are extracted and 95 % confidence level values are marked with dots on the raster plot. This analysis can be performed with a much longer time series, other variables, ad different spatial extent across the globe due to the availability of the GLDAS variables.

Geographic Information System (GIS) analyses are an essential part of natural resource management and research. Calculating and summarizing data within intersecting GIS layers is common practice for analysts and researchers. However, the various tools and steps required to complete this process are slow and tedious, requiring many tools iterating over hundreds, or even thousands of datasets. USGS scientists will combine a series of ArcGIS geoprocessing capabilities with custom scripts to create tools that will calculate, summarize, and organize large amounts of data that can span many temporal and spatial scales with minimal user input. The tools work with polygons, lines, points, and rasters to calculate relevant summary data and combine them into a single output table that can be easily incorporated into statistical analyses. These tools are useful for anyone interested in using an automated script to quickly compile summary information within all areas of interest in a GIS dataset.

Toolbox Use
License
Creative Commons-PDDC
Recommended Citation
Welty JL, Jeffries MI, Arkle RS, Pilliod DS, Kemp SK. 2021. GIS Clipping and Summarization Toolbox: U.S. Geological Survey Software Release. https://doi.org/10.5066/P99X8558

The downloadable ZIP file contains model documentation and contact information for the model creator. For more information, or a copy of the project report which provides greater model detail, please contact Ryan Urie - traigo12@gmail.com.This model was created from February through April 2010 as a central component of the developer's master's project in Bioregional Planning and Community Design at the University of Idaho to provide a tool for identifying appropriate locations for various land uses based on a variety of user-defined social, economic, ecological, and other criteria. It was developed using the Land-Use Conflict Identification Strategy developed by Carr and Zwick (2007). The purpose of this model is to allow users to identify suitable locations within a user-defined extent for any land use based on any number of social, economic, ecological, or other criteria the user chooses. The model as it is currently composed was designed to identify highly suitable locations for new residential, commercial, and industrial development in Kootenai County, Idaho using criteria, evaluations, and weightings chosen by the model's developer. After criteria were chosen, one or more data layers were gathered for each criterion from public sources. These layers were processed to result in a 60m-resolution raster showing the suitability of each criterion across the county. These criteria were ultimately combined with a weighting sum to result in an overall development suitability raster. The model is intended to serve only as an example of how a GIS-based land-use suitability analysis can be conceptualized and implemented using ArcGIS ModelBuilder, and under no circumstances should the model's outputs be applied to real-world decisions or activities. The model was designed to be extremely flexible so that later users may determine their own land-use suitability, suitability criteria, evaluation rationale, and criteria weights. As this was the first project of its kind completed by the model developer, no guarantees are made as to the quality of the model or the absence of errorsThis model has a hierarchical structure in which some forty individual land-use suitability criteria are combined by weighted summation into several land-use goals which are again combined by weighted summation to yield a final land-use suitability layer. As such, any inconsistencies or errors anywhere in the model tend to reveal themselves in the final output and the model is in a sense self-testing. For example, each individual criterion is presented as a raster with values from 1-9 in a defined spatial extent. Inconsistencies at any point in the model will reveal themselves in the final output in the form of an extent different from that desired, missing values, or values outside the 1-9 range.This model was created using the ArcGIS ModelBuilder function of ArcGIS 9.3. It was based heavily on the recommendations found in the text "Smart land-use analysis: the LUCIS model." The goal of the model is to determine the suitability of a chosen land-use at each point across a chosen area using the raster data format. In this case, the suitability for Development was evaluated across the area of Kootenai County, Idaho, though this is primarily for illustrative purposes. The basic process captured by the model is as follows: 1. Choose a land use suitability goal. 2. Select the goals and criteria that define this goal and get spatial data for each. 3. Use the gathered data to evaluate the quality of each criterion across the landscape, resulting in a raster with values from 1-9. 4. Apply weights to each criterion to indicate its relative contribution to the suitability goal. 5. Combine the weighted criteria to calculate and display the suitability of this land use at each point across the landscape. An individual model was first built for each of some forty individual criteria. Once these functioned successfully, individual criteria were combined with a weighted summation to yield one of three land-use goals (in this case, Residential, Commercial, or Industrial). A final model was then constructed to combined these three goals into a final suitability output. In addition, two conditional elements were placed on this final output (one to give already-developed areas a very high suitability score for development [a "9"] and a second to give permanently conserved areas and other undevelopable lands a very low suitability score for development [a "1"]). Because this model was meant to serve primarily as an illustration of how to do land-use suitability analysis, the criteria, evaluation rationales, and weightings were chosen by the modeler for expediency; however, a land-use analysis meant to guide real-world actions and decisions would need to rely far more heavily on a variety of scientific and stakeholder input.

This data package includes spatial environmental and social layers for Shivamogga District, Karnataka, India that were considered as potential predictors of patterns in human cases of Kyasanur Forest Disease (KFD). KFD is a fatal tick-borne viral haemorrhagic disease of humans, that is spreading across degraded forest ecosystems in India. The layers encompass a range of fifteen metrics of topography, land use and land use change, livestock and human population density and public health resources for Shivamogga District across 1km and 2km study grids. These spatial proxies for risk factors for KFD that had been jointly identified between cross-sectoral stakeholders and researchers through a co-production approach. Shivamogga District is the District longest affected by KFD in south India. The layers are distributed as 1km and 2km GeoTiffs in Albers equal area conic projection. For KFD, spatial models incorporating these layers identified characteristics of forest-plantation landscapes at higher risk for human KFD. These layers will be useful for modelling spatial patterns in other environmentally sensitive infectious diseases and biodiversity within the district.

Methods Processing of environmental predictors of Kyasanur Forest Disease distribution

This file details the sources and processing of environmental predictors offered to the statistical analysis in the paper. All processing was performed in the raster package [1] of the R program [2] unless otherwise specified, with function names as specified below.

Topography predictors

Elevation data was extracted in tiles from Shuttle Radar Topography Mission data version 4 [3] an original resolution of 0.000833 degrees Latitude and Longitude resolution (approximately 90m by 90m grid cells). Tiles were mosaicked across the study region using the merge function. A slope value for each pixel was calculated (in degrees) using the terrain function of the raster package, and a focal window of 3 by 3 cells. Both the resulting elevation and slope rasters were cropped to the administrative boundaries of the Shivamogga District (raster package: crop function) and re-projected to an equal area projection (Albers equal area conic projection) using the projectRaster function (method=”bilinear”). Mean elevation and slope values were then calculated across the study 1km and 2km grid cells, using the aggregate function to average values across the appropriate number of ~90m grid cells and then the resample function to align the resulting grid to the study grids.

Landscape predictors

Metrics of the current availability (and fragmentation) of forest, agricultural and built-up land use types as well as that of water-bodies were extracted from the MonkeyFeverRisk Land Use Land Cover map of Shimoga. The latter was produced from classification of earth observation data from 2016 to 2018 using the methods described in the Supplementary information S3 file of the paper linked to this dataset. The LULC map had an original grid square resolution of 0.000269 degrees Latitude and Longitude resolution (or 30m x 28m grid cells) and nine different LULC classes. It was cropped to the administrative boundaries of the Shimoga District (raster package: crop function) and re-projected to the equal area projection (Albers equal area conic projection) using the projectRaster function (method=”ngb” for categorical data). The agriculture and fallow land classes were combined before landscape analysis (due to the difficulty of separating them accurately in the classification process).

An algorithm was developed in R to identify which of the pixels in the LULC map coincided with each 1km and 2km grid cell of the study area. The ClassStat function of the SDM Tools package [4] was used to calculate the proportional area of each 1km or 2km grid cell landscape that was made up of a particular land class, as well patch density and edge density metrics for the forest classes as indicators of fragmentation and forest-agriculture interface habitat respectively (Fig. S2B). The proportional area values (pi) of the n different forest classes (wet evergreen forest, moist deciduous forest, dry deciduous forest and plantation) were used to calculate an index of forest type diversity per grid cell as follows, after Shannon & Weaver (1949) [5]:

H'= -1npi(lognpi)

Metrics of longer term forest changes in Shimoga since 2000 were derived from a global product by Hansen et al. (2013) [6] available at a spatial resolution of 1 arc-second per pixel, (~ 30 meters per pixel at equator). Forest loss during the period 2000–2014, is defined as a stand-replacement disturbance, or a change from a forest to non-forest state, encoded as either 1 (loss) or 0 (no loss). Forest gain during the period 2000–2012, is defined as a non-forest to forest change entirely within the study period, encoded as either 1 (gain) or 0 (no gain).These layers were again cropped to the administrative boundaries of the Shimoga District (raster package: crop function) and re-projected to an equal area projection (Albers equal area conic projection) using the projectRaster function (method=”ngb”) in R. An algorithm was developed in R to identify which of the pixels in the loss and gain rasters coincided with each 1km and 2km grid cell of the study area. The ClassStat function of the SDM Tools package [4] was used to calculate the proportional area of each 1km or 2km grid cell that was made up of loss pixels or gain pixels. Forest gain and loss are very highly correlated (r=0.986) and occur in similar places in the landscape (Fig. S2C). Forest loss was a much more common transition than a forest gain affecting 1.2% of land pixels rather than 0.16% of land pixels for forest gain.

To assess how forest loss or gain from a global product like Hansen et al. (2013) should be interpreted locally in south India, we analysed how the loss and gain pixels from Hansen et al. 2013 coincided with classes in the MonkeyFeverRisk LULC map (by extracting the value of the LULC map for the centroids of loss or gain pixels).

The distribution of loss and gain pixels across forest classes from the MonkeyFeverRisk LULC map is shown in Table 1. Locations categorised as a loss by Hansen et al. were most commonly classified currently as plantation, followed by moist evergreen forest, followed by

moist or dry deciduous forest by the MonkeyFeverRisk LULC map. The pattern was similar for the gain pixels. Since not all forest loss pixels were non-forest in the current day and not all forest gain pixels were forest in the current day, the precise meaning of the Hansen et al. (2013) forest loss layer was unclear for south India, though we expect that it is at least indicative of areas where the forest has undergone a larger degree of change since 2000.

Table 1: Percentage of loss (n= 108398) and gain (n= 14646) land pixels from the global Hansen et al. (2013) product that fall into different forest classes according to the MonkeyFeverRisk LULC map

    Land use class


    Gain


    Loss




    moist evergreen


    30.4


    26.1




    moist deciduous


    6.5


    16.2




    dry deciduous


    3.0


    9.7




    plantation


    46.2


    37.2




    Non-forest classes


    14.0


    10.9

Host and public health predictors

Livestock host density data, namely buffalo and indigenous cattle densities in units of total head per village were obtained from Department of Animal Husbandry, Dairying and Fisheries, Government of India Census from 2011 at village level. These were linked to village boundaries from the Survey of India using the village census codes in R. The village areas were calculated from the spatial polygons dataframe of villages using the rgeos package in R, so that the total head per village metrics could be convert into an areal density of buffalo and indigenous cattle per km and then rasterized at 1km and 2km using the rasterize function of the raster package.

The human population size and public health metrics were obtained from the Government of India Population Census 2011. The human population size (census field TOT_P) was again linked to the spatial polygon village boundaries using the census village code (census field VCT_2011) and converted to an areal metric of population density per km and rasterized at 1km and 2km as above. The number of medics per head of population was derived by summing all doctors and para-medicals “in position” across all types of health centres, clinics and dispensaries per village and dividing by the total population of the village (TOT_P) and then linked to village boundaries and rasterized as above. The proximity to health centres was a categorical variable derived from the “Primary.Health.Centre..Numbers” field, where 1 = Primary Health Centre (PHC) within village boundary, 2 = PHC within 5km of village, 3=PHC within 5-10km of village, 4= PHC further than 10km from village. It was linked to village boundaries and rasterized as above.

The resulting raster layers for all predictors were saved in GeoTiff format.

References

Robert J. Hijmans (2017). raster: Geographic Data Analysis and Modeling. R package version 2.6-7. https://CRAN.R-project.org/package=raster
R Core Team (2018). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria.URL https://www.R-project.org/ 
Jarvis, A., Reuter, I., Nelson, A., Guevara, E. Hole-filled SRTM for the globe Version 4. 2008.
VanDerWal, J., Falconi, L., Januchowski, S., Shoo, L., & Storlie, C. (2014). SDMTools: Species Distribution Modelling Tools: Tools for processing data associated with species distribution modelling exercises. R

scripts.zip

arcgisTools.atbx: terrainDerivatives: make terrain derivatives from digital terrain model (Band 1 = TPI (50 m radius circle), Band 2 = square root of slope, Band 3 = TPI (annulus), Band 4 = hillshade, Band 5 = multidirectional hillshades, Band 6 = slopeshade). rasterizeFeatures: convert vector polygons to raster masks (1 = feature, 0 = background).

makeChips.R: R function to break terrain derivatives and chips into image chips of a defined size. makeTerrainDerivatives.R: R function to generated 6-band terrain derivatives from digital terrain data (same as ArcGIS Pro tool). merge_logs.R: R script to merge training logs into a single file. predictToExtents.ipynb: Python notebook to use trained model to predict to new data. trainExperiments.ipynb: Python notebook used to train semantic segmentation models using PyTorch and the Segmentation Models package. assessmentExperiments.ipynb: Python code to generate assessment metrics using PyTorch and the torchmetrics library. graphs_results.R: R code to make graphs with ggplot2 to summarize results. makeChipsList.R: R code to generate lists of chips in a directory. makeMasks.R: R function to make raster masks from vector data (same as rasterizeFeatures ArcGIS Pro tool).

terraceDL.zip

dems: LiDAR DTM data partitioned into training, testing, and validation datasets based on HUC8 watershed boundaries. Original DTM data were provided by the Iowa BMP mapping project: https://www.gis.iastate.edu/BMPs. extents: extents of the training, testing, and validation areas as defined by HUC 8 watershed boundaries. vectors: vector features representing agricultural terraces and partitioned into separate training, testing, and validation datasets. Original digitized features were provided by the Iowa BMP Mapping Project: https://www.gis.iastate.edu/BMPs.

This collection of files are part of a larger dataset uploaded in support of Low Temperature Geothermal Play Fairway Analysis for the Appalachian Basin (GPFA-AB, DOE Project DE-EE0006726). Phase 1 of the GPFA-AB project identified potential Geothermal Play Fairways within the Appalachian basin of Pennsylvania, West Virginia and New York. This was accomplished through analysis of 4 key criteria: thermal quality, natural reservoir productivity, risk of seismicity, and heat utilization. Each of these analyses represent a distinct project task, with the fifth task encompassing combination of the 4 risks factors. Supporting data for all five tasks has been uploaded into the Geothermal Data Repository node of the National Geothermal Data System (NGDS).

This submission comprises the data for Thermal Quality Analysis (project task 1) and includes all of the necessary shapefiles, rasters, datasets, code, and references to code repositories that were used to create the thermal resource and risk factor maps as part of the GPFA-AB project. The identified Geothermal Play Fairways are also provided with the larger dataset. Figures (.png) are provided as examples of the shapefiles and rasters. The regional standardized 1 square km grid used in the project is also provided as points (cell centers), polygons, and as a raster. Two ArcGIS toolboxes are available: 1) RegionalGridModels.tbx for creating resource and risk factor maps on the standardized grid, and 2) ThermalRiskFactorModels.tbx for use in making the thermal resource maps and cross sections. These toolboxes contain item description documentation for each model within the toolbox, and for the toolbox itself. This submission also contains three R scripts: 1) AddNewSeisFields.R to add seismic risk data to attribute tables of seismic risk, 2) StratifiedKrigingInterpolation.R for the interpolations used in the thermal resource analysis, and 3) LeaveOneOutCrossValidation.R for the cross validations used in the thermal interpolations.

Some file descriptions make reference to various 'memos'. These are contained within the final report submitted October 16, 2015.

Each zipped file in the submission contains an 'about' document describing the full Thermal Quality Analysis content available, along with key sources, authors, citation, use guidelines, and assumptions, with the specific file(s) contained within the .zip file highlighted.

UPDATE: Newer version of the Thermal Quality Analysis has been added here: https://gdr.openei.org/submissions/879 (Also linked below) Newer version of the Combined Risk Factor Analysis has been added here: https://gdr.openei.org/submissions/880 (Also linked below) This is one of sixteen associated .zip files relating to thermal resource interpolation results within the Thermal Quality Analysis task of the Low Temperature Geothermal Play Fairway Analysis for the Appalachian Basin. This file contains 6 images (.png) including predicted and associated error for surface heat flow, depth to 80 degrees C, depth to 100 degrees C, temperature at 1.5 km, temperature at 2.5 km and temperature at 3.5 km.

The sixteen files contain the results of the thermal resource interpolation as binary grid (raster) files, images (.png) of the rasters, and toolbox of ArcGIS Models used. Note that raster files ending in “pred” are the predicted mean for that resource, and files ending in “err” are the standard error of the predicted mean for that resource. Leave one out cross validation results are provided for each thermal resource.

Several models were built in order to process the well database with outliers removed. ArcGIS toolbox ThermalRiskFactorModels contains the ArcGIS processing tools used. First, the WellClipsToWormSections model was used to clip the wells to the worm sections (interpolation regions). Then, the 1 square km gridded regions (see series of 14 Worm Based Interpolation Boundaries .zip files) along with the wells in those regions were loaded into R using the rgdal package. Then, a stratified kriging algorithm implemented in the R gstat package was used to create rasters of the predicted mean and the standard error of the predicted mean. The code used to make these rasters is called StratifiedKrigingInterpolation.R Details about the interpolation, and exploratory data analysis on the well data is provided in 9_GPFA-AB_InterpolationThermalFieldEstimation.pdf (Smith, 2015), contained within the final report.

The output rasters from R are brought into ArcGIS for further spatial processing. First, the BufferedRasterToClippedRaster tool is used to clip the interpolations back to the Worm Sections. Then, the Mosaic tool in ArcGIS is used to merge all predicted mean rasters into a single raster, and all error rasters into a single raster for each thermal resource.

A leave one out cross validation was performed on each of the thermal resources. The code used to implement the cross validation is provided in the R script LeaveOneOutCrossValidation.R. The results of the cross validation are given for each thermal resource.

Other tools provided in this toolbox are useful for creating cross sections of the thermal resource. ExtractThermalPropertiesToCrossSection model extracts the predicted mean and the standard error of predicted mean to the attribute table of a line of cross section. The AddExtraInfoToCrossSection model is then used to add any other desired information, such as state and county boundaries, to the cross section attribute table. These two functions can be combined as a single function, as provided by the CrossSectionExtraction model.

This geodatabase includes spatial datasets that represent the Southeastern Coastal Plain aquifer system in the States of Alabama, Georgia, Mississippi, South Carolina, and Tennessee. Included are: (1) polygon extents; datasets that represent the aquifer system extent, the entire extent subdivided into subareas or subunits, and any polygon extents of special interest (outcrop areas, no data available, areas underlying other aquifers, anomalies, for example), (2) raster datasets for the altitude of each aquifer subarea or subunit, (3) altitude, and/or if applicable, thickness contours used to generate the surface rasters, (4) georeferenced images of the figures that were digitized to create the altitude and thickness contours. The images and digitized contours are supplied for reference. The extent of the Southeastern Coastal Plain aquifer system is derived the linework in the Southeastern Coastal Plain aquifer system extent maps in a digital version of the aquifer extent presented in the Groundwater Atlas of the United States (the U.S. Geological Survey Hydrologic Atlas HA-730-F, -730-G, and -730-K. The Southeastern Coastal Plain aquifer system has 4 aquifer subunits, in order from the most surficial to the deepest: A1: Chickasawhay River aquifer, A2: Pearl River Aquifer, A3: Chattahoochee River Aquifer, and A4: Black Warrior River Aquifer. The altitude and thickness contours for each available subunit were digitized from georeferenced figures of altitude contours in U.S. Geological Survey Profession Paper 1410-B, (USGS PP 1410-B), and the resultant top and bottom altitude values were interpolated into surface rasters within a GIS using tools that create hydrologically correct surfaces from contour data, derive the altitude from the thickness (depth from the land surface), and merge the subareas into a single surface. The primary tool was "Topo to Raster" used in ArcGIS, ArcMap, Esri 2014. The surface rasters were corrected for the areas where the altitude of an underlying layer of the aquifer exceeded altitude of an overlying layer.

Here we present three datasets describing three large European landscapes in France (Bauges Geopark - 89,000 ha), Poland (Milicz forest district - 21,000 ha) and Slovenia (Snežnik forest - 4,700 ha) down to the tree level. Individual trees were generated combining inventory plot data, vegetation maps and Airborne Laser Scanning (ALS) data. Together, these landscapes (hereafter virtual landscapes) cover more than 100,000 ha including about 64,000 ha of forest and consist of more than 42 million trees of 51 different species. For each virtual landscape we provide a table (in .csv format) with the following columns:- cellID25: the unique ID of each 25x25 m² cell- sp: species latin names- n: number of trees. n is an integer >= 1, meaning that a specific set of species "sp", diameter "dbh" and height "h" can be present multiple times in a cell.- dbh: tree diameter at breast height (cm)- h: tree height (m) We also provide, for each virtual landscape, a raster (in .asc format) with the cell IDs (cellID25) which makes data spatialisation possible. The coordinate reference systems are EPSG: 2154 for the Bauges, EPSG: 2180 for Milicz, and EPSG: 3912 for Sneznik. The v2.0.0 presents the algorithm in its final state. Finally, we provide a proof of how our algorithm makes it possible to reach the total BA and the BA proportion of broadleaf trees provided by the ALS mapping using the alpha correction coefficient and how it maintains the Dg ratios observed on the field plots between the different species (see algorithm presented in the associated Open Research Europe article). Below is an example of R code that opens the datasets and creates a tree density map. ------------------------------------------------------------# load package library(terra) library(dplyr)

set work directory

setwd() # define path to the I-MAESTRO_data folder

load tree data

tree <- read.csv2('./sneznik/sneznik_trees.csv', sep = ',')

load spatial data

cellID <- rast('./sneznik/sneznik_cellID25.asc')

set coordinate reference system

Bauges:

crs(cellID) <- "epsg:2154"

Milicz:

crs(cellID) <- "epsg:2180"

Sneznik:

crs(cellID) <- "epsg:3912"

convert raster into dataframe

cellIDdf <- as.data.frame(cellID) colnames(cellIDdf) <- 'cellID25'

calculate tree density from tree dataframe

dens <- tree %>% group_by(cellID25) %>% summarise(n = sum(n))

merge the two dataframes

dens <- left_join(cellIDdf, dens, join_by(cellID25))

add density to raster

cellID$dens <- dens$n

plot density map

plot(cellID$dens)

Link to the ScienceBase Item Summary page for the item described by this metadata record. Service Protocol: Link to the ScienceBase Item Summary page for the item described by this metadata record. Application Profile: Web Browser. Link Function: information

This raster maps average atmospheric methane concentrations in 2016. The data used to create it is from the NASA Aqua satellite (specifically the AIRS instrument) that records monthly average atmospheric methane concentrations. AIRS collects methane data at different pressure levels. The raster depicts data at the 400 hPa level because that is where the instrument is most sensitive to methane concentration. The monthly data was consolidated using the NASA tool, Giovanni https://giovanni.gsfc.nasa.gov/giovanni/, to create a raster with annual average methane concentrations. Giovanni output two rasters: one for daytime averaged data and one for nighttime averaged data. ArcGIS was then used to combine the two rasters to create a single annual raster. A conversion factor of 1.0e+9 was multiplied to convert the final raster from mole fractions to parts per billion. The results are attached. The Carnegie Endowment for International Peace would like to eventually incorporate a methane raster as a new layer in the Oil Climate Index (OCI) web tool http://oci.carnegieendowment.org/. Carnegie is currently updating the OCI, adding greenhouse gas comparisons of global gas fields and visualizing their methane emissions. Carnegie is planning to work with our OCI partners at Stanford to further analyze the methane concentration raster to separate out signal from noise and to identify potential methane concentration hot spots associated with oil and gas operations. This raster is useful when studying short term climate risks, especially when it comes to Arctic oil and gas resources.

This hosted feature layer has been published in RI State Plane Feet NAD 83. Conservation ecologists have coined the term Ecological Land Units (ELU) to describe and map the physical properties of landscapes. Typically, ELUs are defined by the geology, soils, elevation, and landform (hilltop, hillside, valley). A specific ELU has a unique combination of soils, geology, landform, and elevation. ELUs are derived from soil and elevation data using a GIS. It was important that we used readily available data and we kept the derivation of ELUs as simple as possible. After consulting the published literature and conferring with expert soil scientists and plant ecologists, we focused on two aspects of soils, soil drainage class and soil texture. Soil drainage class is very good at distinguishing wet versus dry habitats. Soil texture (sandy, silty, loamy, etc.) is an important habitat component for plants. Using USDA SSURGO (State Soil Survey Geographic Database) data that is readily available from RIGIS, we created a raster dataset (50 feet cell size) of the different soil drainage classes and another raster dataset of the soil texture classes. There are many properties of soils that are available to use for analyses such as this, for example stoniness, depth to bedrock, etc. The two factors we chose are extremely important soil properties in supporting different plant communities. Landform represents where a location is with respect to elevation, slope, and aspect (direction a hillside is facing). Landform distinguishes hilltops, hill sides, valley bottoms, etc. We used the RIGIS digital terrain model as our source of elevation data to measure landform. Landform classes were identified using GIS modeling of slope, aspect, and elevation. The final ELU map is made by adding together the raster datasets for landform, drainage class, and soil texture. Because we were careful with our encoding system, the sum of the three rasters provides us a composite of the individual datasets. For example, a location that is a well-drained (code value 2000) and consists of gravelly sand (code value 100) a sits on a hilltop (code value 21) and would combine to be ELU 2121 (2000+100+21). This process yielded 204 unique ELUs for the state of Rhode Island. Examination of a cumulative distribution function (CDF) of the ELUs showed that most of the ELUs were small and did not occur very often. Conversely, 20 ELUs were quite common and encompassed almost 85% of the land area of RI.Find out more about Mapping ELUs

Alberta Satellite Land Cover (ASLC) raster is produced for Government of Alberta and is available to the general public. The ASLC Mosaic is a land cover dataset for the province of Alberta. It is a set of individual projects completed in the last decade for wildfire historical burns and various areas. As the remote sensing technology and related information extraction approaches have made significant advancements during last 20 years the process naturally evolved from classification schema targeting Landsat 5 to one integrating various sources of remote sensing data, information extraction procedures, different open and commercial toolboxes for processing and more automation. Integration of multiple data types includes various optical and radar data free of charge as well as GOA's purchased LiDAR data swaths. Since raster values in products from individual projects are identical with values of initial Alberta Ground Cover Classification product, it was possible to merge all those into new ASLC raster.

These contiguous suitable area data are based on the land suitability data from the Southern Gulf Water Resource Assessment in the NT and Qld (SOGWRA). To address operational farming constraints imposed by parcels of suitable land being too small or oddly shaped according to natural variability of land, or physical limits on suitable farming land parcel sizes, contiguous suitable area data was generated. This contiguous suitable area data is based on crop-specific minimum areas and minimum length/width of contiguous suitable land and is produced as standalone data products for all crop groups. The rules are provided for download. The data was generated to remove the component of landscape complexity that natural distributions of soil and land variability and specific crop requirements produce. This data provides improved land evaluation information to identify opportunities and promote detailed investigation for a range of sustainable development options. The land suitability evaluation methods used to produce the underlying data are a modification of the Food and Agriculture Organisation (FAO) land evaluation approach. The land suitability analysis is described in full in the CSIRO SOGWRA published report ‘Soils and land suitability for the Southern Gulf catchments’. A technical report from the CSIRO Southern Gulf Water Resource Assessment to the Government of Australia. The naming convention for these data is; ‘crop group’ underscore ‘major crop’ underscore ‘season code’ underscore ‘irrigation type code’ underscore ‘catchments code’ underscore ‘data type’ eg ‘CG7_CottonGrains_D_Fw_S_ContigArea’ is Cotton and grain crops grown in the dry season with furrow irrigation in the Southern Gulf catchments contiguous suitable area data. The codes for season are; W – wet season; D – dry season; P – perennial. The codes for irrigation type are; S – overhead spray irrigation; T – trickle irrigation; Fd – flood irrigation; Fw – furrow irrigation; R – rainfed. It is important to emphasize that this is a regional-scale assessment: further data collection and detailed soil physical, chemical and nutrient analyses would be required to plan development at a scheme, enterprise or property scale. Several limitations that may have a bearing on land suitability were out of scope and not assessed as part of this activity (refer to the report), these limitations include biophysical and socio-cultural. For example these land suitability raster datasets do not include consideration of the licensing of water, flood risk, contiguous land, risk of irrigation induced secondary salinity, or land tenure and other legislative controls. Some of these may be addressed elsewhere in SOGWRA eg flooding was investigated by the Earth observation remote sensing group in the surface water activity. The Southern Gulf Water Resource Assessment provides a comprehensive overview and integrated evaluation of the feasibility of aquaculture and agriculture development in the Southern Gulf catchments (NT and Qld) as well as the ecological, social and cultural (indigenous water values, rights and aspirations) impacts of development. Lineage: These contiguous suitable area raster datasets have been generated from a range of inputs and processing steps. Following is an overview. For more information refer to the CSIRO SOGWRA published report ' Soils and land suitability for the Southern Gulf catchments’. A technical report from the CSIRO Southern Gulf Water Resource Assessment to the Government of Australia. 1. Collated existing data (relating to: soils, climate, topography, natural resources, remotely sensed, of various formats: reports, spatial vector, spatial raster etc). 2. Selection of additional soil and land attribute site data locations by a conditioned Latin hypercube statistical sampling method applied across the covariate data space. 3. Fieldwork was carried out to collect new attribute data, soil samples for analysis and build an understanding of geomorphology and landscape processes. 4. Database analysis was performed to extract the data to specific selection criteria required for the attribute to be modelled. 5. The R statistical programming environment was used for the attribute computing. Models were built from selected input data and covariate data using predictive learning from a Random Forest approach implemented in the ranger R package. 6. Create Digital Soil Mapping (DSM) attribute raster datasets. DSM data is a geo-referenced dataset, generated from field observations and laboratory data, coupled with environmental covariate data through quantitative relationships. It applies pedometrics - the use of mathematical and statistical models that combine information from soil observations with information contained in correlated environmental variables, remote sensing images and some geophysical measurements. 7. Land management options were chosen and suitability rules created for DSM attributes. 8. Suitability rules were run to produce limitation subclass datasets using a modification on the FAO methods. 9. Final suitability data created for all land management options. 10. Companion predicted reliability data was produced from the 500 individual Random Forest attribute models created. 11. QA Quality assessment of these land suitability data was conducted by two methods. Method 1: Statistical (quantitative) assessment of the "reliability" of the spatial output data presented as a raster of the Confusion Index. Method 2: A workshop was conducted in March 2023 to review DSM soil attribute and land suitability products and facilitated an alternative to the field external validation carried out in other northern Australia water resource assessments. Stakeholders from the NT and Qld jurisdictions reviewed, evaluated and discussed the soundness of the data and processes. The workshop desk top assessment approach provided recommendations for acceptance, improvement and re-modelling of attributes based on expert knowledge and understanding of the soil distribution and landscape in the study area and available data. 12. A two-step process was developed to simplify the data and was applied across the suitability data of the catchments. First the five suitability classes were aggregated to two: ‘suitable’ for suitability classes 1, 2 and 3, or ‘not suitable’ for class 4 and 5. Second, to further simplify the data, and to reflect the on-ground spatial constraints of farming practices, isolated one or two pixels of ‘not suitable’ contained in larger ‘suitable’ areas were reclassified as ‘suitable’. 13. For each crop group, a minimum area and width were defined based on knowledge of farming practices. Depending on the possible land use, minimum areas were deemed as 2.5 ha, 5 ha, 10 ha or 25 ha and minimum widths of 80 m or 120 m (rules are provided for download). 14. For each crop rule the minimum width was imposed by removing those parts of the suitable area that are narrower (in any direction) than the required minimum width. The remaining groups of connected cells were then tested to see if they meet the required minimum area and removed if they did not.

This map is the subset of the World Terrestrial Ecosystems map, prepared specifcally for the Pacific Region. The World Terrestrial Ecosystems map classifies the world into areas of similar climate, landform, and land cover, which form the basic components of any terrestrial ecosystem structure. This map is important because it uses objectively derived and globally consistent data to characterize the ecosystems at a much finer spatial resolution (250-m) than existing ecoregionalizations, and a much finer thematic resolution (431 classes) than existing global land cover products.Cell Size: 250-meter Source Type: ThematicPixel Type: 16 Bit UnsignedData Projection: GCS WGS84Extent: GlobalSource: USGS, The Nature Conservancy, EsriUpdate Cycle: NoneWhat can you do with this layer?This map allows you to query the land surface pixels and returns the values of all the input parameters (landform type, landcover/vegetation type, climate region) and the name of the terrestrial ecosystem at that location.This layer can be used in analysis at global and local regions. However, for large scale spatial analysis, we have also provided an ArcGIS Pro Package that contains the original raster data with multiple table attributes. For simple mapping applications, there is also a raster tile layer. This layer can be combined with the World Protected Areas Database to assess the types of ecosystems that are protected, and progress towards meeting conservation goals. The WDPA layer updates monthly from the United Nations Environment Programme.Developing the World Terrestrial EcosystemsWorld Terrestrial Ecosystems map was produced by adopting and modifying the Intergovernmental Panel on Climate Change (IPCC) approach on the definition of Terrestrial Ecosystems and development of standardized global climate regions using the values of environmental moisture regime and temperature regime. We then combined the values of Global Climate Regions, Landforms and matrix-forming vegetation assemblage or land use, using the ArcGIS Combine tool (Spatial Analyst) to produce World Ecosystems Dataset. This combination resulted of 431 World Ecosystems classes.Each combination was assigned a color using an algorithm that blended traditional color schemes for each of the three components. Every pixel in this map is symbolized by a combination of values for each of these fields.The work from this collaboration is documented in the publication:Sayre et al. 2020. An assessment of the representation of ecosystems in global protected areas using new maps of World Climate Regions and World Ecosystems - Global Ecology and Conservation More information about World Terrestrial Ecosystems can be found in this Story Map.

Reason for SelectionIslands provide important habitat for many species, including birds, sea turtles, mammals, insects, and plants. Their relative isolation from disturbance and mainland predators can make them important breeding habitat for coastal birds and sea turtles (as represented by piping plover and loggerhead sea turtle). Their unique ecology and isolation can also make them important habitat for some mammals, plants and insects that are only found on islands (as represented by Cape Sable thoroughwort, Florida semaphore cactus, silver rice rat, and Bartram’s hairstreak butterfly). As a barrier that can protect the mainland from major storms, they also help protect ecosystems and human communities from extreme weather events. The critical habitat included in this indicator refers to areas with specific physical or biological features that are essential to conserving a federally threatened or endangered species and may require special management or protection.Input DataBase Blueprint 2022 extentSoutheast Blueprint 2023 extentIsland boundaries from the Global Island Explorer provided by the U.S. Geological Survey and Esri, accessed 5-13-2022. (Global Islands data has been updated since this release. Click here to download the latest version.)Critical habitat provided by the U.S. Fish and Wildlife Survey, accessed 6-23-2022Mapping StepsClip the Global Island Explorer data to the Base Blueprint 2022 extent and merge small and large islands.From the critical habitat data (CRITHAB_POLY.shp), select the following species based on ‘comname’: piping plover, loggerhead sea turtle, Cape Sable thoroughwort, Florida semaphore cactus, silver rice rat, Bartram’s hairstreak butterfly.Clip the critical habitat for the selected species to the merged islands layer.Convert the islands and selected species data to raster and clip to the spatial extent of Base Blueprint 2022.Add zero values to help users better understand the extent of this indicator and to make it perform better in online tools. Buffer the island shapefile by 40 km to make a continuous buffer along the coast, with no gaps. Use this to create a raster of zeros for that buffer.Combine rasters so parts of islands with selected species critical habitat get a value of 2, other island pixels get a value of 1, and all other areas in the buffer described above get a value of 0.As a final step, clip to the spatial extent of Southeast Blueprint 2023.Note: For more details on the mapping steps, code used to create this layer is available in the Southeast Blueprint Data Download under > 6_Code. Final indicator valuesIndicator values are assigned as follows:2 = Island critical habitat for any of six threatened and endangered species (piping plover, loggerhead sea turtle, Cape Sable thoroughwort, Florida semaphore cactus, silver rice rat, or Bartram’s hairstreak butterfly)1 = Other island area0 = Not a coastal islandKnown IssuesThis indicator underestimates piping plover critical habitat in parts of the Chandeleur Islands off the coast of Louisiana. Locations of barrier islands can be highly dynamic, and the island boundaries and critical habitat data did not agree on the locations of some parts of the islands.Disclaimer: Comparing with Older Indicator Versions There are numerous problems with using Southeast Blueprint indicators for change analysis. Please consult Blueprint staff if you would like to do this (email hilary_morris@fws.gov).Literature Cited Sayre, R., S. Noble, S. Hamann, R. Smith, D. Wright, S. Breyer, K. Butler, K. Van Graafeiland, C. Frye, D. Karagulle, D. Hopkins, D. Stephens, K. Kelly, Z, basher, D. Burton, J. Cress, K. Atkins, D. van Sistine, B. Friesen, B. Allee, T. Allen, P. Aniello, I Asaad, M. Costello, K. Goodin, P. Harris, M. Kavanaugh, H. Lillis, E. Manca, F. Muller-Karger, B. Nyberg, R. Parsons, J. Saarinen, J. Steiner, and A. Reed. 2018. A new 30 meter resolution global shoreline vector and associated global islands database for the development of standardized global ecological coastal units. Journal of Operational Oceanography–A Special Blue Planet Edition. [https://doi.org/10.1080/1755876X.2018.1529714]. U.S. Geological Survey and Esri. Global Island Explorer. Accessed May 13, 2022. [https://rmgsc.cr.usgs.gov/gie/]. U.S. Fish and Wildlife Service. Critical Habitat. Accessed June 23, 2022. [https://ecos.fws.gov/ecp/report/table/critical-habitat.html].

This land suitability for Banana raster data (in GeoTIFF format) represents areas of potential suitability for this crop and its specific irrigation management systems in the Flinders and Gilbert catchments of North Queensland. The data is coded 1-5: 1 - Suitable with no limitations; 2 - Suitable with minor limitations; 3 - Suitable with moderate limitations; 4 - Marginal; 5 - Unsuitable. The land suitability evaluation methods used to produce this data are a modification of methods of the Food and Agriculture Organisation of the UN (FAO). This data is part of the Flinders and Gilbert Agricultural Resource Assessment (FGARA) project and is designed to support sustainable regional development in Australia being of importance to Australian Governments and agricultural industries. The project identifies new opportunities for irrigation development in these remote areas by providing improved soil and land evaluation data to identify opportunities and promote detailed investigation. A companion dataset exists, “Confidence of suitability data for the FGARA project”. A link to this dataset can be found in the “related materials” section of this metadata record. Lineage: These suitability raster data for Banana and its individual irrigation management systems have been created from a range of inputs and processing steps. Below is an overview. For more information refer to the CSIRO FGARA published reports and in particular: Bartley R, Thomas MF, Clifford D, Phillip S, Brough D, Harms D, Willis R, Gregory L, Glover M, Moodie K, Sugars M, Eyre L, Smith DJ, Hicks W and Petheram C (2013) Land suitability: technical methods. A technical report to the Australian Government for the Flinders and Gilbert Agricultural Resource Assessment (FGARA) project, CSIRO. Broadly, the steps were to: 1. Collate existing data (data related to: climate, topography, soils, natural resources, remotely sensed etc of various formats; reports, spatial vector, spatial raster etc). 2. Select additional soil and attribute site data by Latin hypercube statistical sampling method applied across the covariate space. 3. Carry out fieldwork to collect additional soil and attribute data and understand geomorphology and landscapes. 4. Build models from selected input data and covariate data using predictive learning via rule ensembles in the RuleFit3 software. 5. Create Digital Soil Mapping (DSM) key attributes output data. DSM is the creation and population of a geo-referenced database, generated using field and laboratory observations, coupled with environmental data through quantitative relationships. It applies pedometrics - the use of mathematical and statistical models that combine information from soil observations with information contained in correlated environmental variables, remote sensing images and some geophysical measurements. 6. Choose land management options and create suitability rules for DSM attributes. 7. Run suitability rules to produce limitation datasets using a modification on the FAO methods. 8. Create final suitability data for all land management options. Two companion datasets exist for this dataset. The first is linked to in the “related materials” section of this metadata record, entitled “Confidence of suitability data for the FGARA project”. The second (held by CSIRO Land and Water) includes expert opinion and knowledge about landscape processes or conditions that will influence agricultural development potential in these catchments, but were not captured sufficiently in the modelling process (and areas of expert opinion where the Mahanabolis method underestimates confidence). The two landscape features that require special attention are the basalt rock outcrops in the Upper Flinders catchment that were not well captured by the covariate data, and the secondary salinisation hazard in the central Flinders catchment. For more information refer to the report “Land suitability: technical methods. A technical report to the Australian Government for the Flinders and Gilbert Agricultural Resource Assessment (FGARA) project”.

The Southeast Blueprint uses a least-cost path connectivity analysis to identify connections between priority areas. A program called Linkage Mapper defines corridors that link hubs across the shortest distance possible, while also routing through as much Blueprint priority as possible. The hubs are defined as large patches of highest priority Blueprint areas and/or protected areas. Because the continental and Caribbean portions of the Southeast Blueprint are non-contiguous, we run a separate connectivity analysis for each area. We use the same overall methods in both places; the primary difference is a smaller size threshold for hubs in the Caribbean (500 acres vs. 5,000 acres) to reflect the relative maximum size potential compared to the continental Southeast. Other minor differences between the continental and Caribbean methods—like the use of wildlife road crossing data where available, removal of reservoirs to align with the Blueprint methods, and running the analysis in a buffer area around the Blueprint extent to align with neighboring plans—are detailed in the subsequent mapping steps for each analysis. The combined results of the connectivity analysis, the Southeast Blueprint 2024 hubs and corridors, are available on the Blueprint page of the SECAS Atlas. CONTINENTAL CONTINENTAL RESISTANCE RASTER & HUBS This is the resistance raster or cost surface used in the Linkage Mapper-based connectivity analysis for the continental portion of Blueprint 2024. Input Data

2024 Southeast Blueprint subregions 2024 Southeast Blueprint combined Zonation results 2023 Midwest Conservation Blueprint combined Zonation results 2023 Midwest Conservation Blueprint potential hubs (provided by Rachael Carlberg with MLI on 7-24-2024; for more information contact rachael_carlberg@fws.gov) The Nature Conservancy’s (TNC) Terrestrial Resilience - Local Connectedness (Anderson et al. 2016) Florida Wildlife Road crossings, last updated 12-7-2020 Select North Carolina wildlife road crossings

Eastern NC crossings were provided by Gary Jordan at the Raleigh, NC U.S. Fish and Wildlife Service Ecological Services field office on 7-22-2021; Gary received these data from Travis Wilson at the North Carolina Wildlife Resources Commission. Other wildlife road crossings in NC were provided by Alex Vanko from the Wildlands Network and Steve Goodman from the National Parks Conservation Association. Some of these are dedicated crossing structures and some are bridge replacements that were specifically widened to accommodate wildlife; this is not a census of all widened bridge replacements.

OpenStreetMap data “roads” layer, accessed 2-23-2022. These were used to estimate the fencing that is built around the wildlife road crossings.

OpenStreetMap® is open data, licensed under the Open Data Commons Open Database License (ODbL) by the OpenStreetMap Foundation (OSMF). Additional credit to OSM contributors. Read more on the OSM copyright page.

2021 National Land Cover Database (NLCD) U.S. Geological Survey (USGS) Watershed Boundary Dataset (WBD), accessed 12-2-2021: HUC2s; download the data Protected Areas Database of the United States (PAD-US): PAD-US 3.0 national geodatabase - Combined Proclamation Marine Fee Designation Easement; PAD-US 4.0 national geodatabase - Combined Proclamation Marine Fee Designation Easement Reservoirs used to remove highly altered areas in inland continental subregions (available in the Southeast Blueprint 2024 Ancillary Data Download)

Mapping Steps Define Extent of the Continental Connectivity Analysis

Dissolve the continental Southeast subregions into a single polygon covering the 2024 continental Blueprint extent. Buffer the continental extent by 26.1 km. The distance of 26.1 km is based on dispersal distances of subadult black bears (White et al. 2000), a species that disperses from within the Southeast Blueprint geography into neighboring areas. Convert the buffered continental extent from vector to raster.

Begin Creation of the Resistance Raster

The resistance raster is primarily derived from three data sources: the Southeast Blueprint combined Zonation results, the Midwest Blueprint potential hubs, and TNC local connectedness. Wildlife road crossings and protected areas will be integrated in a later step. Flip the Southeast Blueprint continental combined Zonation results so that areas that ranked higher in Zonation are easier to move through. We already converted the Zonation from a floating point raster to an integer raster in a previous step while rebalancing to account for reservoirs that were not included in the prioritization.To match the Southeast results, convert the Midwest Blueprint combined Zonation results to integer and flip so that areas that ranked higher in Zonation are easier to move through. These values are used in the buffer area around the continental Southeast Blueprint to allow corridors to connect into neighboring landscapes. Rescale and flip the local connectedness layer from TNC’s Resilient Land project so it has similar values to the Zonation results. Clip these data to the buffer area around the Southeast Blueprint. Using the ArcPy Rescale “Linear” function, rescale the local connectedness layer to assign values ranging from 1 to 100. Also, flip the high and low values so that areas with high local connectedness are easier to move through (i.e., lower resistance) than areas of low local connectedness. Convert these results to integer. This will be used in the buffer area around the continental Southeast Blueprint in places that aren't covered by the Midwest Blueprint. Mosaic together the resulting Zonation outputs and the TNC local connectedness layer, using the Southeast Zonation results anywhere they occur, the Midwest Zonation results in places where there are no Southeast Zonation results, and the local connectedness in the remaining buffer area.

Make Urban Areas Harder to Move Through

To get the Linkage Mapper connectivity analysis to run successfully across the whole area, we had to resample the resistance raster from 30 m to 90 m. During this process, developed areas often get resampled out. To make sure that the connectivity analysis has sufficient information to move around developed areas, select the high and medium developed classes from the 2021 NLCD, resample them to 90 m, and burn them into the resistance raster with a high resistance value (150).

Make Known Wildlife Road Crossings Easier to Move Through We obtained a small sample of known wildlife road crossings for North Carolina and Florida. This served as a test to route corridors through these crossings. We would like to get more spatial data for wildlife road crossings in the future.

Buffer known wildlife road crossings by 210 m and convert them to raster with a low resistance value (1) that makes it easy for corridors to route through. Select roads surrounding known wildlife road crossings from the OpenStreetMap data. We do not want corridors to go across these roads except at the designated crossings, so we needed to make the resistance very high in these areas. Buffer the roads by 180 m and convert to raster with a very high resistance value (150). In testing, the methods above weren’t enough to encourage corridors to route through the crossings, so we buffered the crossings again, this time by 900 m, erased the road areas from the larger buffered crossings (so the roads remained a high value) and converted to a raster with value that is easier to move through (10). Combine the roads with the combined Zonation results, using a maximum value so that these roads are very hard to cross. Combine the above raster with the two buffered road crossing layers, using a minimum value so that the areas with wildlife road crossings get the lower value from the buffered road crossing layers or Zonation results.

Make Protected Areas Easier to Move Through

Protected conservation lands are an important component of regional connectivity. Many are captured as hubs in a later step, but not all conserved areas meet the size threshold to qualify as a hub. To better capture their value, regardless of size, burn protected areas into the resistance raster to make them easier to move through. To do this, first prepare the protected areas data, starting with the PAD-US 4.0 combined proclamation, marine, fee, designation, and easement layer. To exclude areas that do not meet the intent of this step to capture the value of protected conservation lands for connectivity:Remove areas with location designations of ‘School Trust Land’, ‘School Lands’, ‘School Land’, ‘State Land Board’, or ‘3201’. These extensive lands are leased out and are not open to the public. Remove areas with the designation type of 'Military' or “Proclamation'. Military lands are not primarily managed for conservation. The proclamation category represents the approved boundary of public lands, within which land protection is authorized to occur, but not all lands within the proclamation boundary are necessarily currently in a conserved status. Remove areas with the owner name of 'BOEM'. These Outer Continental Shelf lease blocks do not represent actively protected marine parks, but serve as the “legal definition for BOEM offshore boundary coordinates...for leasing and administrative purposes” (BOEM). Remove areas with a category of 'Proclamation' (see explanation above).

PAD-US 4.0 is missing state wildlife management area boundaries in Oklahoma. Extract those from PAD-US 3.0 by using a combination of a state name of 'Oklahoma' and local designation of 'State Wildlife Management Area'. Merge the selected polygons from PAD-US 4.0 and PAD-US 3.0, then convert to raster and reclassify to assign a value of 0 to unprotected areas and a value of 10 to protected areas. Subtract the protected area raster from the resistance raster. Lowering the resistance value makes it easier for corridors to route

Reason for Selection Beaches in Puerto Rico and the U.S. Virgin Islands support a wide array of shorebirds, colonial seabirds, and sea turtles (ACJV et al. 2015, USFWS 2022). However, their limited spatial extent makes beaches an ecosystem of special concern for conservation (ACJV et al. 2015). In addition, coastal dunes and beaches are some of the Caribbean ecosystems facing the greatest threat from disturbance and development (ACJV et al. 2015). This indicator focuses on a suite of bird and sea turtle species that nest on beaches, though it includes beach habitat used for other activities like foraging and breeding, in addition to nesting. Input Data

  Puerto Rico Gap Analysis Project predicted vertebrate species distributions: data provided by Dr. Bill Gould with the Caribbean Climate Hub on 4-4-2022 (contact william.a.gould@usda.gov for more information); read the final report 

In Puerto Rico, we used the following GAP species models: Wilson’s plover

Puerto Rico Gap Analysis Project landcover; download the data; read the final report 
U.S. Virgin Islands Gap Analysis Project predicted vertebrate species distributions and landcover; data and report appendices provided by Dr. Bill Gould with the Caribbean Climate Hub on 2-6-2023 (contact william.a.gould@usda.gov for more information); read the final report 

In the USVI, we used the following GAP species models:Wilson’s ploverAmerican oystercatcherHawksbill sea turtleLeatherback sea turtleGreen sea turtle

State of the World’s Sea Turtles (SWOT) nest locations for hawksbill, leatherback, green, and loggerhead sea turtles; download the data using the download icon on the left side of the SWOT mapping application. Note: loggerhead was only observed in USVI, while the other species were observed in both PR and USVI. 
OpenStreetMap data, accessed 6-28-2023 

A polygon from this dataset is considered a beach if the value in the “nature” attribute is beach. OpenStreetMap describes natural areas as “a wide variety of physical geography, geological and landcover features”. Data were downloaded in .pbf format using R code. OpenStreetMap® is open data, licensed under the Open Data Commons Open Database License (ODbL) by the OpenStreetMap Foundation (OSMF). Additional credit to OSM contributors. Read more on the OSM copyright page.

U.S. Virgin Islands beaches layer (separate vector layers for St. Croix, St. Thomas, and St. John) provided by Joe Dwyer with Lynker/the NOAA Caribbean Climate Adaptation Program on 3-3-2023 (contact jdwyer@lynker.com for more information)
Southeast Blueprint 2023 subregions: Caribbean 
Southeast Blueprint 2023 extent

Mapping Steps

Buffer by 15 m the beach polygons from OpenStreetMap and the USVI beaches layer. This is consistent with the methodology used in the urban park size indicator to avoid the loss of narrow beaches when converting to raster. Project and convert to raster.
Extract from GAP landcover the relevant beach classes and resample to a 30 m pixel. In Puerto Rico, the beaches include two classes: “Gravel beaches and stony shoreline” and “Fine to coarse sandy beaches, mixed sand and gravel beaches”. In the U.S. Virgin Islands, the beaches include three classes: “Fine to Medium Grained Sandy Beaches,” “Gravel Beaches” and “Mixed Sand and Gravel Beaches”. 
Extract the predicted habitat class from the GAP predicted habitat rasters for species that nest on beaches. These include Wilson’s plover, American oystercatcher, and hawksbill, leatherback, and green sea turtles. Note: Only Wilson’s plover was predicted in Puerto Rico. Project and resample the rasters to 30 m. 
Extract from the SWOT data nest all point locations for hawksbill, leatherback, green and loggerhead sea turtles and convert them to 30 m rasters. 
Merge together the SWOT and GAP predicted habitat rasters for each species and identify each pixel that contains at least one species. Then clip the resulting raster to the beach extent.
To define individual beaches, run a region group on beach extent.
Run the Zonal Statistics “MAX” function to apply species presence to the entirety of each beach.
Reclassify to 0 the beach extent layer created above.
To create the final indicator values seen below, mosaic together three rasters: the beaches containing at least one species, the beach extent, and the Caribbean Blueprint 2023 extent. This adds back in a 0 value for areas not identified as beaches and a 1 value for beaches that did not contain any species predictions or observations. Zero values better represent the extent of the source data and make the indicator perform better in online tools.
Clip to the Caribbean Blueprint 2023 subregion.
As a final step, clip to the spatial extent of Southeast Blueprint 2023. 

Note: For more details on the mapping steps, code used to create this layer is available in the Southeast Blueprint Data Download under > 6_Code. Final indicator values Indicator values are assigned as follows: 2 = Beach with 1+ nesting species predicted or observed 1 = Other beach 0 = Not identified as a beach Known Issues

This indicator likely underprioritizes beaches in Puerto Rico due to disparities in both beach extent and species data coverage between Puerto Rico and the U.S. Virgin Islands. USVI has a comprehensive hand-digitized beach layer that is not available in Puerto Rico. GAP models only one beach-nesting species in Puerto Rico, compared to five in the U.S. Virgin Islands (though the additional species are known to occur in Puerto Rico as well). In addition, the SWOT sea turtle observations better aligned with the beach polygons in USVI. We will explore additional datasets and methods for addressing these disparities in future revisions. To help mitigate this issue for this year, we set the maximum species richness as 1+ rather than using the full range of species richness values, since Puerto Rico had a maximum species richness of 3 in the available data, compared to a maximum value of 6 in USVI.
This indicator includes beaches from OpenStreetMap, which is a crowdsourced dataset. While members of the OpenStreetMap community often verify map features to check for accuracy and completeness, there is the potential for spatial errors (e.g., misrepresenting the boundary of a beach) or incorrect tags (e.g., labelling an area as a beach that is not actually a beach). However, using a crowdsourced dataset gives on-the-ground experts, Blueprint users, and community members the power to fix errors and add new parks to improve the accuracy and coverage of this indicator in the future.
This indicator may exclude some small beaches that aren’t captured in the source data. We encourage interested partners and citizens to add any missing beaches to OpenStreetMap so we can better capture them in future updates.
This indicator does not account for other factors that influence the quality of beach habitat, such as distance to roads, light pollution, and vulnerability to erosion and sea-level rise.

Other Things to Keep in Mind

The species chosen for this indicator are birds and sea turtles that nest on beaches. However, the indicator also includes beach habitat used for activities other than nesting, like foraging and breeding.
This indicator does not always align with Caribbean coastal shoreline condition. Some areas identified as important beach habitat in this indicator, especially those coming from the GAP Wilson’s plover model, are scored as armored in coastal shoreline condition (e.g., the Hyatt Regency Grand Hotel in Río Grande, Puerto Rico). This often occurs where riprap is present along narrow beaches, or occasionally near bulkheads. There is often a section of beach present behind the riprap or bulkhead that could still provide habitat, or the riprap is sporadically placed on a long stretch of beach to protect inland structures. In these cases, the mismatch reflects the different intent of these complementary indicators. In some cases, hardened structures may be actually misclassified as beach. Inconsistencies in alignment and classification likely result from the older age and coarser resolution of the GAP data (10 m raster based on 2001 landcover) compared to the more recent and fine-scale CUSP shorelines (vectors dating primarily from 2014-2021) and challenges in distinguishing the unique remote sensing signature of beach vs. riprap and other hardened structures. Because of the 30 m resolution of the Blueprint and underlying data, a single pixel may contain a mix of beach habitat and hardened structures and be reflected differently in each of these two indicators due to their different functions.

Disclaimer: Comparing with Older Indicator Versions There are numerous problems with using Southeast Blueprint indicators for change analysis. Please consult Blueprint staff if you would like to do this (email hilary_morris@fws.gov). Literature Cited Atlantic Coast Joint Venture, Caribbean Landscape Conservation Cooperative, and U.S. Fish and Wildlife Service. Avian Conservation Planning Priorities for Puerto Rico and the U.S. Virgin Islands (BCR 69). February 2015. [https://acjv.org/documents/PRUSVI_plan.pdf].

Gould, William A.; Alarcón, Caryl; Fevold, Brick; Jiménez, Michael E.; Martinuzzi, Sebastián; Potts, Gary; Quiñones, Maya; Solórzano, Mariano; Ventosa, Eduardo. 2008. The Puerto Rico Gap Analysis Project. Volume 1: Land cover, vertebrate species distributions, and land stewardship. Gen. Tech. Rep. IITF-GTR-39. Río Piedras, PR: U.S. Department of Agriculture, Forest Service, International Institute of Tropical Forestry. 165 p. [https://data.fs.usda.gov/research/pubs/iitf/iitf_gtr39.pdf].

Gould WA, Solórzano MC, Potts GS, Quiñones M, Castro-Prieto J, Yntema LD. 2013. U.S. Virgin Islands Gap Analysis Project – Final Report. USGS, Moscow ID and the USDA FS International Institute of