Study site and population. Koram Island (12.242°, 100.009°) lies ~ 1 km offshore in the Gulf of Thailand and within Khao Sam Roi Yot National Park, Prachuap Khiri Khan, Thailand. It has an area of 0.45 km2 and a coastline of 3.5 km. The habitat––limestone karst blanketed with a dense flora of dwarf evergreen trees and deciduous scrub, and encircled by rocky shore and sandy beaches––supports a population of ca. 75 long-tailed macaques described as hybrids at the subspecies taxonomic level (Macaca fascicularis aurea × M. f. fascicularis)30. The animals are well habituated to human observers as a result of regular tourism and sustained study since 201331.
Rank determination. Macaques form multi-male multi-female social groups with individual dominance hierarchies. Among females, this hierarchy is strictly linear and stable through time32. To determine the rank-order of adults, we recorded dyadic agonistic interactions and their outcomes (i.e., aggression, supplants, and silent-bared-teeth displays of submission) during 5-min focal follows of individuals based on a randomised order of continuous rotation31. In some cases, these data were supplemented with ad libitum observations. This protocol existed during five years (2013–2018) of continual observations before we conducted our experiment in July-August 2018. To determine the effects of dominance rank on individual food-cleaning propensities, and to standardise the data by sex, we followed the methods of Levy et al.33. We calculated ordinal ranks between 1 (highest) and n (lowest), where n is the number of animals aged ≥ 5 years in each hierarchy. We analysed male (n = 8) and female (n = 16) dominance hierarchies together in the same statistical models. We chose ordinal ranks because they best approximate competitive regimes during density-dependent competition33. This type of feeding competition reflects those of our experimental trials and the natural conditions on Koram Island, where preferred resources are limited34.
Measuring sand. To quantify the amount of sand on food surfaces, whether provisioned by tourists (cucumbers, melon, pineapple) or used in our experiment (sliced cucumbers), we applied a quick-drying liquid polymer––granulated plastic (Pioloform BL 16; Wacker-Chemie GMBH, Munich, Germany) mixed with ethanol (18% plastic to 82% ethanol, by weight)––to each food item. When dried, we peeled and stored the sand-infused film for analysis. The advantages of this method are twofold: the removal of exogenous particulate matter is extremely thorough; and the plastic does not detach biogenic silica such as trichromes or phytoliths35. In the lab, we dissolved each peel in ethanol and separated the sand by centrifugation, producing a pellet. We dried and weighed the pellet, dividing the mass by the surface area of the food object, which we calculated from digital photographs imported into ImageJ v. 1.52. This method produces an estimate of exogenous particulate mass per area (mg mm− 2), allowing direct comparison of apples and oranges.
To measure the elemental and physical properties of sand, we dispersed and filtered the pellets in water using a 0.2-µm isopore membrane filter, which we submitted for scanning electron microscopy (SEM) and energy dispersive X-ray spectroscopy (EDS). To establish the parameters for multi-field and bulk analysis, we imaged a representative area of the filter at multiple magnifications and performed discrete particle analysis. 50X magnification allowed for statistically significant representation of particle number and size range (allowing a 5 µm lower particle size range in analysis). All discrete particle analyses indicated silicon rich particles and composition distribution bins were established to include dominant accompanying elements. After establishing these parameters, we initiated multi-field automated analysis using six fields of view of the debris field (at 50X). A composition classification was assigned to each particle and data sorted by composition classification and particle size (particle sizing was binned using standard Feret maximum parameter). The sizing bins are standard ISO-16232 size classes.
Experimental design. To elicit food handling behaviours and determine individual cleaning preferences, we put three cucumber slices in each of three trays (20 x 30 x 10 cm) and manipulated the amount of contaminating sand. In the low-sand treatment, we put cucumber slices in a tray without sand; however, the presence of some aeolian sand was unavoidable. In the intermediate-sand treatment, we lined the tray with sand and put the cucumber slices on the surface. In the high-sand treatment, we buried cucumber slices completely (Fig. 2a). Trays were placed 1.5 m apart and 15 m from the ocean, and we randomised the colour and sequence of trays across trials.
Trials began when one or more monkeys approached the trays and ended when the animals finished every cucumber slice or abandoned the experiment (range: 10 s to 14 min). We used video recordings to determine the onset and offset of individual food-handling bouts, beginning from initial contact with a cucumber slice and ending when the final slice, in its entirety, entered the mouth. Within each bout, we determined the duration of brushing and washing behaviours, defining each from the onset of serial stereotypical forelimb movements to the moment of oral ingestion. We estimated energy intake rates by calculating the number of cucumber slices consumed during each food-handling bout, and multiplying each slice by 1.1 kcal (source: U.S. Department of Agriculture, FoodData Central, 2019). We performed 101 trials over 5 weeks and recorded 1,282 food-handling bouts by 42 individual monkeys. To minimise the potential confounding effects of dominance interactions, we analysed trials with ≤ 3 monkeys. Thus, 935 food-handling bouts were included in GLMM statistical models, which included data on individual rank, sex, and sand treatment. If a monkey consumed a cucumber slice without brushing or washing it, the zero-second duration was included in both GLMMs.
Behavioural analyses. To model experimental variance in brushing and washing behaviours as a function of experimental treatment, sex, and rank, we fit generalised linear mixed models (GLMMs) using the glmmTMB package in R (version 4.2.3)36. In each GLMM, we modelled sand treatment (a categorical variable with three levels), sex (a categorical variable with two levels), and ordinal rank (a discrete variable ranging from 1 to 18) as fixed effects. We also incorporated two additional interaction terms as fixed effects: sand treatment × ordinal rank and sand treatment × sex. To account for experimental variance among individuals and control for pseudoreplication37 (because the number of feeding bouts per individual varied widely), we included individual ID as a random intercept. The brushing and washing datasets were whole-number counts (seconds) with means < 5. The distributions were right-skewed with high concentrations of biologically-meaningful zeros38 (i.e., instances of food-handling without any cleaning behaviour). Thus, we fit four separate models in glmmTMB to account for these non-normal spreads: standard GLMMs with a log-link function and either a (1) Poisson or (2) negative binomial error distribution (default = nbinom2 in glmmTMB); or zero-inflated generalised linear mixed models (ZIGLMM) with a logit-link function, a single zero-inflation parameter applying to all observations, and either a (3) Poisson or (4) negative binomial error distribution. We then determined the best fit model using delta AIC values in the bbmle package in R.
dAIC values indicated that, for both brushing (Table S1) and washing (Table S2), negative binomial GLMMs without zero-inflation were the best-fitting models. We validated each model by calculating dispersion statistics (𝜒2/degrees of freedom)39. Dispersion statistics for the brushing model (𝜒2 = 541.92; residual degrees of freedom = 564; 𝜒2/rdf = 0.96, p = 0.74, one-sided test) and the washing model (𝜒2 = 174.22; residual degrees of freedom = 351; 𝜒2/rdf = 0.50, p = 1.00, one-sided test) failed to detect overdispersion in either case. We report the fixed effects tests for each GLMM in Tables S3 and S4 as Analysis of Deviance Tables (Type II Wald chi square tests, one sided) along with 𝜒2 values, degrees of freedom, and p-values (one-sided tests). For all statistical analyses, α = 0.05.
Optimal cleaning time model. To model the optimisation of sand removal before consumption, we accounted for two distinct temporal periods: handling time \(h\), which includes an assessment time and pre-cleaning time, and the cleaning time \(t\). Assessment time (set as a constant 1 second) includes visual fixation on a food object and forelimb extension before contact, whereas pre-cleaning time represents all handling activities that precede cleaning. During brushing, the pre-cleaning time was essentially nil (zero seconds), but washing required travel from the experimental treatments to the ocean, requiring longer pre-cleaning times (x̄ = 22 ± 15 s; range: 5 to 78 s). We assumed that the proportion of sand removed from each cucumber follows the saturating relationship \(g\left(t\right) = t/(c+t)\), where \(c\) is the half-saturation constant associated with brushing or washing. As \(c\) increases, so does the inefficiency of a given cleaning behaviour. Given our observations that brushing removes 75% of grit in 2.97 s, and washing removes 93% of grit in 3.53 s (Fig. S1), we obtain the constants \({c}_{brushing}= 0.99\) s and \({c}_{washing}=0.26\) s, such that washing (without considering handling costs) is the most efficient strategy. The rate of grit removal is then given by \(R\left(t\right) = g\left(t\right)/(h + t)\), which reaches a maximum at the optimal cleaning time \({t}^{*}= \sqrt{c h}\). For brushing and washing cleaning strategies, we obtain the expected optimal cleaning times \({{t}^{*}}_{brushing}=0.98\) s, and \({{t}^{*}}_{washing}=2.39\) s (Fig. 3a), respectively. These optimal cleaning times are defined exclusively with respect to maximising the rate of grit removal, without considering the potentially cascading effects of these strategies on fitness.