South Arabia’s prehistoric monument landscape shows social resilience to climate change

PLOS One | https://doi.org/10.1371/journal.pone.0323544 May 28, 2025 7 / 24 initial construction our studies of stone placements, arrangements, inclusions, additions, and radiocarbon dates indicate. As we documented through excavation, the reuse of HCT monuments occurred sometimes millennia after construction and disuse. This circumstance extended the HCT timeline in Hadramawt (pink dotted line in Fig 2). Elsewhere in Oman, Williams and Gregoricka [40] document re-visits and deposits subsequent to SCAB construction. Our SCAB sample (n = 23) is undated, but we find some that share attributes of Neolithic platforms, such as a standing stone and East- facing façade (e.g., D106-002); others are only some meters from Neolithic platforms (e.g., D027-001, D027-003). (S2 Table, S3 Table) Environmental data analysis and machine learning To determine the factors influencing monument placement, we performed bootstrap aggregating of classification trees and regularized multinomial logistic regression. We extracted cell values from 12 geospatial raster datasets at the 371 monu- ment locations and from a random selection of 500 true absence locations within the viewshed. The environmental raster layers are listed in Table 1 and hosted in the PANGAEA data repository [78]. Topographic variables were derived from the elevation layer using the Geomorphometry and Gradient Metrics toolbox [64] in ArcGIS [79]. We calculated path distances between monuments and each of human routes and springs. Path distances (anisotro-

pic) accounted for cost distances in travelling upslope versus downslope. The cost raster generated in GIS, and used to calculate path distances, used a symmetric inverse linear as the vertical factor setting, which sets an exponentially increasing cost with slope steepness. We identified springs through place name lists and intensive field coverage targeting Table 1. The independent environmental variables used in the analysis. Independent variable Description of variable Eastness Aspect decomposed to two directions, with highest values on slopes facing east and lowest values on slopes facing west. Elevation An elevation raster layer derived from ASTER Global Digital Elevation Map V2. The original data has been reprojected and resampled using bilinear interpolation. Fog density A multiband (R, G, B, NIR) raster layer of spatial variability in monsoon fog density calculated on a per cell basis as the mean of the fog reflectance values of 119 Landsat 5 TM scenes, 17 Landsat 7 ETM+ scenes and 121 Landsat 8 OLI TIRS scenes. Normalized Difference Vegetation Index (NDVI) Maximum Greenest Pixel A max-NDVI raster layer where pixel value is the maximum greenest pixel of all NDVI layers derived from Landsat 8. Northness Aspect decomposed to two directions, with highest values on slopes facing north and lowest values on slopes facing south. Path distance to roads Path distances were calculated from a cost raster in ArcGIS. A cost raster identifies the cost of travelling through each cell. In our cost raster, the cost of travel increased exponentially with slope steepness.

Path distance to springs Slope [75] A raster layer of slope, in degrees, derived from the elevation layer using the average maximum technique. Terrain Ruggedness Index (TRI) [76] A TRI raster layer derived from the elevation layer, based on the sum change in elevation between a central cell and its eight neighboring cells. TRI is a terrain roughness metric. Topographic Position Index (TPI) (R = 200, 1000, 5000 m) (Guisan et al., 1999) [77] TPI raster layers with radii (R) of 200 m, 1 km and 5 km, derived from the elevation layer. TPI measures the difference between a central cell eleva- tion and the average elevation around it within a predetermined radius (R). At small scales it is a measure of terrain roughness and at large scales it describes slope position and landform types. https://doi.org/10.1371/journal.pone.0323544.t001