Multiproxy approach to reconstruct fossil primate feeding behavior: Case study for macaque from the Plio-Pleistocene site Guefaït-4.2 (eastern Morocco)

The genus Macaca belongs to Cercopithecidae (Old World monkeys), Cercopithecinae, Papionini. The presence of Macaca in North Africa is well known from the Late Miocene to the Late Pleistocene. However, the diet of fossil Macaca has been poorly described in the literature. In this study, we investigated the feeding habits of Macaca cf. sylvanus (n = 4) from the Plio-Pleistocene site Guefaït-4.2 in eastern Morocco through multiproxy analysis combining analyses of stable carbon and oxygen isotopes from tooth enamel, buccal microtexture, and low-magnification occlusal dental microwear. For both microwear analyses, we compared the macaques with a new reference collection of extant members of Cercopithecoidea. Our occlusal microwear results show for the fossil macaque a pattern similar to the extant Cercocebus atys and Lophocebus albigena, African forest-dwelling species that are characterized by a durophagous diet based mainly on hard fruit and seed intake. Buccal microtexture results also suggest the consumption of some grasses and the exploitation of more open habitats, similar to that observed in Theropithecus gelada. The δ13C of M. cf. sylvanus indicates a C3 based-diet without the presence of C4 plants typical of the savanna grassland in eastern Africa during this period. The high δ18O values of M. cf. sylvanus, compared with the contemporary ungulates recovered from Guefaït-4.2, could be associated with the consumption of a different resource by the primate such as leaves or fresh fruits from the upper part of trees. The complementarity of these methods allows for a dietary reconstruction covering a large part of the individual’s life.


Introduction
Modern Cercopithecoidea primates inhabit a wide range of habitats, including tropical forests, high-altitude woodlands, and savanna grasslands (Jablonski et al., 2000). This adaptability across a range of environments began with a radiation within and beyond Africa, probably due to the increasingly open and fragmented environments of the Late Miocene and the Plio-Pleistocene (Elton, 2007;Hughes et al., 2008;Roos et al., 2019). For example, the closing of the Tethys Ocean and subsequent aridification of northern Africa (Zhang et al., 2014), and fluctuations in Milankovitch cycles (DeMenocal, 1995), formed intermittent biogeographical barriers. It is hypothesized that these ecological barriers diversified Cercopithecoidea, specifically the Papionini tribe, into different genera. Theropithecus adapted to wet lowlands, Papio, Cercocebus, Lophocebus, and Mandrillus filled niches within more forested environments, while Macaca became more-generalized and settled a wide variety of habitats within higher latitudes (Delson, 1980). Paleontological and molecular evidence indicates that the macaque lineage originated during the Miocene in Africa and diverged from other African Papionini around 7-5.5 Ma (Delson, 1980;Tosi et al., 2005;Roos et al., 2019). According to the fossil record, macaques dispersed throughout Eurasia during the Late Miocene, inhabiting primarily the Mediterranean region 6-5 Ma (Alba et al., 2014;Piñero et al., 2017) and then in the Late Pliocene, eastern Asia (~4 Ma; Szalay and Delson, 1979;Elton, 2007;Meloro and Elton, 2012). In this radiation, macaques experienced a rapid diversification into several species (Elton and O'Regan, 2014;Roos et al., 2019). The genus Macaca consists of approximately 23 extant species and occupies more diverse ecological environments and geographical ranges (from 10° south to over 40° north latitude) than any other genus of nonhuman primates. In fact, the expansion of this genus in Africa and Eurasia is one of the most successful among primates (Jablonski, 2002;Jablonski and Frost, 2010;Fleagle, 2013;Roos et al., 2019), and its ecological distribution and adaptative capacity reflect great dietary diversity and interspecific behavioral responses to seasonality (Kato et al., 2014).
At present, there is only one species of the genus Macaca in Africa, M. sylvanus (Barbary macaques) from Algeria (Chiffa, Tigounatine, Icetcifère, Akfadou, and Kherrata) and Morocco (High Atlas, Middle Atlas, and Rif;Fooden, 2007). Field observations have defined feeding habits as varied and to include an important seasonal component, with diet consisting of more than 100 species of plants, fungi, lichens, mosses, and animal prey (Deag, 1983;Ménard, 1985;Ménard and Vallet, 1986). Studies that have taken into account the annual composition of the diet indicate that seeds, leaves, and herbs are the main foods of Barbary macaques in Akfadou and Tigounatine (Algeria; Ménard, 1985;Ménard andVallet, 1986, 1996). However, in the Rif mountains of northern Morocco, more than 50% of their diet is based on seeds (Mehlman, 1988), and some studies of the populations that inhabit the Middle Atlas suggest that fruits could constitute a significant part of their diet (Deag, 1983;Drucker, 1984;Ménard and Mohamed, 1999), demonstrating high intra-species local variation at the population level. Nevertheless, the diet of fossil M. cf. sylvanus is virtually unknown, with the exception of a study on dental microwear of Macaca sp. from the Pliocene of northwestern China, which suggests a hard-object feeder or extractive forager, without grass consumption (Williams and Holmes, 2011).
Many cercopithecoid fossils have been recovered from the same sites where Hominin fossils have been found (Delson, 1980). The genus Homo experienced a similar out-of-Africa dispersion during the Pleistocene and inhabited the same ecosystems as members of Papionini. Therefore, both lineages had to confront the same climatic constrictions and habitat partitioning. Understanding the ecology of the fossil Cercopithecoidea primates is of interest for interpreting the paleoecology and niche partitioning of both groups. Consequently, we aimed to analyze the feeding strategies of M. cf. sylvanus from GFT-4.2 on different timescales using, for the first time, a multi-proxy approach. On the one hand, the isotopic signal (δ 18 O and δ 13 C) of tooth enamel bioapatite provides dietary information during the first years of the animal's life following the tooth mineralization times proposed by Sirianni and Swindler (1985) for the extant Macaca nemestrina. On the other hand, dental microwear analysis, combining buccal (microtexture analysis) and occlusal (low magnification) surfaces, provides dietary information over the long term (weeks to Frontiers in Ecology and Evolution 03 frontiersin.org months; Romero et al., 2012;Hernando et al., 2021) and short term (last days; Semprebon et al., 2004), respectively, but always toward the end of the animal's life. Occlusal microwear studies in primates have suggested that mechanical properties of food items are the primary cause of different microwear pattern among primates. For instance, primates that incorporate hard objects, such as seeds and nuts, have more pitted surfaces than folivores primates that exhibit heavily scratched areas (Ungar, 1996;Scott et al., 2006). Occlusal microwear has a short-term turnover therefore, occlusal microwear patterns provide information on the nature of food intake during the last days prior to death (Teaford and Oyen, 1989). On the other hand, buccal surfaces should only interact with the food items consumed while they are being processed in the oral cavity (Ungar and Teaford, 1996;Romero et al., 2012), and experimental studies have shown that the buccal pattern is more stable through time, at least in the absence of significant dietary shift (Romero et al., 2012). We hypothesize that combining occlusal and buccal microwear analyses will strengthen the dietary discrimination among these species. Based on the dietary information of extant M. sylvanus from Morocco and Algeria, we expect, for the GFT-4.2 specimens, a buccal and occlusal pattern similar to folivore-frugivore primates with the incorporation of seeds during the dry season. In the case of stable carbon and oxygen isotopes, we expect values similar to those of the Theropithecus from Ahl al Oughlam (Morocco). That is, a predominance of C 3 plants in the diet and high oxygen values in relation to ungulates. To compare the dental microwear of M. cf. sylvanus from GFT-4.2, we present a new reference collection of wild modern cercopithecines including the same individuals for both microwear methods. The integration of the data obtained from these three high-resolution techniques has allowed us to obtain a more precise and robust picture of the feeding habits of M. cf. sylvanus from GFT-4.2.
Site:  An interdisciplinary Spanish-Moroccan project started in 2006 as an international cooperation between the Université Mohamed Premier Oujda (UMP, Eastern Morocco) and the Catalan Institute for Human Palaeoecology and Social Evolution from Tarragona (IPHES-CERCA, Catalonia, Spain). The aim of this cooperation is to establish the geochronological and archeo-paleontological context of the Aïn Beni Mathar/Guefaït fluvio-lacustrine basin, in the northern part of the High Plateaus (Jerada Province, Eastern Morocco; Aouraghe et al., 2016;Chacón et al., 2016;Sala-Ramos et al., 2022). One of the systematic surveys in the region revealed the fossiliferous level of Guefaït-4. This horizon yielded a rich and diverse faunal assemblage of vertebrates, which was systematically excavated at the GFT-4.2 site (reptiles, amphibians, and both small and large mammals), including the macaque fossil remains analyzed in this work (Alba et al., 2021). The fossil remains were concentrated in a palustrine level constituted of clay and marls, at the base of Unit 2, within the Dhar Iroumyane stratigraphic section (Aouraghe et al., 2019b). The fieldwork carried out between 2017 and 2019 recovered more than 3,000 fossil remains in an excavated area of 28 m 2 . Biostratigraphic correlations with other sites in Northern Africa, such as the lower levels of the Ain Boucherit site in Algeria (Sahnouni et al., 2018) and Ahl al Oughlam in Morocco (Geraads, 2006), suggest an age close to the Plio-Pleistocene boundary (2.58 Ma). The absence of the genus Equus (despite the abundance of Hipparion equid remains; Aouraghe et al., 2019a,b) indicates an age older than ~2.4 Ma (Sahnouni et al., 2018). This chronological limit is supported by the study of the murid assemblage Piñero et al., 2019), the carnivore assemblage (Madurell-Malapeira et al., 2021), and preliminary paleomagnetism studies (Parés et al., 2020).

Materials
A total of six teeth of M. cf. sylvanus were recovered from GFT-4.2 during the fieldwork in 2018 and 2019 (Alba et al., 2021). The materials are housed in the Faculty of Science, Mohammed I University in Oujda, Morocco. We discarded two of them for the final analysis because they correspond to dental germ fragments (GFT4.2018-1-P13-191, right M2 or M3, and GFT4.2019-1-R15-45, right M2).
The studied teeth (n = 4) correspond to a left upper central incisor (I 1 , GFT4.2018-1-Q14-70), a right lower first molar (M 1, the lingual crown fragment of GFT4.2018-1-Q14-13), a right second lower molar (M 2 , GFT4.2019-1-S15-65), and a left lower third molar (M 3 , GFT4.2019-1-R13-63). All four available fossil teeth (M 1 , M 2 , M 3 , and I 1 ) were sampled for stable isotope analysis, three (M 1 -GFT140, M 2 -GFT138, and M 3 -GFT137) were sampled for buccal dental microtexture analysis, and two (M 1 -GFT140 and M 2 -GFT138) were samples for occlusal dental microwear analysis. M 3 was discarded for the occlusal dental microwear analysis due to the absence of wear facets on this surface. Metrical comparisons suggests that the incisor would belong to a different individual than the lower molars. However, it cannot be ruled out that the lower molars were from the same individual (Alba et al., 2021).
The extant Cercopithecoidea comparative sample for the buccal microtexture and occlusal microwear analyses correspond to a total of 41 well-preserved left second lower molars (M 2 ; Supplementary Table 1). We excluded 10 samples from occlusal microwear analysis and 12 samples from buccal microtexture due to bad preservation or other taphonomic defects on the analyzed surface. The selected species have well-defined habitat and dietary preferences and include five extant species (Theropithecus gelada, Mandrillus sphinx, Lophocebus albigena, Cercocebus atys, and Macaca mulatta). Whenever possible, we chose specimens of the same taxa from close localities to minimize the variability of their dietary regime. We selected the extant species based on the criteria of phylogenetic affinity (Perelman et al., 2011) and the dental morphology that characterizes all cercopithecines, namely bilophodont molars (Kay and Hylander, 1978;Monson and Hlusko, 2014). This comparative sample served to establish a broad spectrum of ecological and dietary diversity.
1. Theropithecus gelada specimens from Ethiopia have a plantbased dietary regime (Dunbar and Bose, 1991), with more than 80% consisting of blades (Mau et al., 2009). However, those T. gelada living on the Guassa Plateau also incorporate large quantities of forbs into their diet (Fashing et al., 2014), and are often selective of tender forbs that are rare on the plateau (Souron, 2018 mechanically challenging items like seeds and hard-shell fruits from the forest floor (Hoshino, 1985;Lahm, 1986;Tutin et al., 1997;Powell et al., 2017). 3. Cercocebus atys forages on the forest's floor and shrub layers, relying on hard-shelled nuts like Sacoglottis gabonensis yearround. Yet, due to seasonal fluctuations, the total contribution of S. gabonensis can range anywhere from 19% to 80% of the monthly diet McGraw et al., 2014). 4. Lophocebus albigena is classified as a hard object feeder (Lambert et al., 2004;McGraw et al., 2011) that consumes a significant portion of hard seeds and nuts (Lambert et al., 2004). In actuality, due to seasonal fluctuations, L. albigena is a generalist that consumes a high (40%) percentage of soft ripe fruit (Lambert et al., 2004) in addition to seeds, flowers and young leaves (Poulsen et al., 2001). During seasons with low fruit availability, L. albigena will eat bark that, although hard, is rich in protein and soluble sugars (Rogers et al., 1994;Olupot, 1998). This fall-back feeding strategy has also been described for Barbary macaques (Macaca sylvanus; Ménard and Mohamed, 1999). 5. Macaca mulatta has the largest geographical distribution with its ranges spanning from Afghanistan to India, Pakistan, China, Thailand and Vietnam, among others (Sengupta et al., 2014). Its wide distribution reflects highly opportunistic feeding strategy, even though fruits may constitute as much as 70% of their diet in some regions (Fooden, 2007), yet less than 9% in other ecosystems, such as Pakistan (Goldstein and Richard, 1989). They are semi-terrestrial primates that can also feed up with species that are in different strata of the forest (Albert et al., 2013).
We have not included extant M. sylvanus specimens because they have a narrow distribution, and some populations are fed by humans like those of Gibraltar (Spain). Therefore, the baseline samples included in the present study consist of wild specimens with well-defined diet and habitat. To contextualize the oxygen and carbon values of M. cf. sylvanus, we have also included contemporaneous herbivore samples from GFT-4.2 corresponding to the same paleontological level 1. These samples include the remains of Tragelaphini (n = 4), Gazella sp. (n = 4), and Hipparion sp. (n = 4). The list of specimens sampled and analyzed and the raw data for stable isotope analysis are provided in Supplementary Table 2.

Stable isotope analysis
Tooth surfaces were cleaned with a tungsten abrasive drill bit to remove any adhering external material. Enamel powder for bulk isotope analysis was then obtained using a diamond-tipped drill, passed along the full length of the buccal surface to ensure a representative measurement for the entire period of enamel formation.
Powdered enamel samples (2-6 mg) were chemically treated at the Biomolecular Laboratory of IPHES. Chemical treatment of samples was based on protocols originally proposed by Koch et al. (1997) and later modified by Tornero et al. (2013). The samples were treated for 4 h in 0.1 M acetic acid [CH 3 COOH] (0.1 mL solution/0.1 mg of sample), neutralized with distilled water, and freeze-dried.
The δ 18 O and δ 13 C values of tooth enamel bioapatite were measured using an automated carbonate preparation device (KIEL-III) coupled to a gas-ratio mass spectrometer (Finnigan MAT 252) at the Environmental Isotope Laboratory, University of Arizona, United States. Samples were reacted with dehydrated phosphoric acid under a vacuum at 70°C. The accuracy and precision of the measurements were checked and calibrated using the calcium carbonate international standards . Replicate measurements of the standards during analysis had errors of ±0.1‰ for δ 18 O and ±0.08‰ for δ 13 C (1 sigma). Carbon isotope composition is reported in δ notation, where δ 13 C = [(Rsample/Rstandard) -1], R = 13 C/ 12 C, and δ is expressed as per mil notation (‰). Carbon and oxygen values are reported relative to PDB (Vienna Pee Dee Belemnite).
In terrestrial ecosystems there is a marked isotopic distinction between the two photosynthetic pathways in plants, C 3 (Calvin cycle) and C 4 (Hatch-Slack cycle), which differ mainly in their discrimination against 13 C during the fixation of CO 2 (Farquhar et al., 1989). Modern plants that use the C 3 photosynthetic pathway (including most trees, woody shrubs, bushes, herbs and temperate grasses) have an average isotopic composition of −28.5‰ (global range −20 to −37‰), with −23‰ as the maximum value recommended for typical C 3 plants and a cutoff of −31.5‰ for closed-canopy forests (Kohn, 2010). Due to the fact that we discarded the presence of that environment in GFT-4.2. For this study, we use the average δ 13 C value of −27‰ for modern C 3 plants, excluding analyses from the understory of closedcanopy forests below −31.5‰ (Kohn, 2010). Modern plants that use the C 4 photosynthetic pathway (primarily tropical grasses and sedges) have an average isotopic composition of −12‰ (global range −10 to −14‰; Cerling and Harris, 1999;Kohn and Cerling, 2002). Crassulacean Acid Metabolism (CAM) plants can have δ 13 C values corresponding to both C 3 and C 4 plants. However, such plants do not seem to be consumed by primates (Cerling et al., 2013a). The δ 13 C values of atmospheric CO 2 have fluctuated over time and this impacts δ 13 C values in plants (Fiedli et al., 1986;Keeling et al., 2017). We have corrected by ~ +1.5 the δ 13 C values of modern plants, due to the difference between the pre-Industrial δ 13 C CO2 value of ~ −6.5‰ established for the Plio-Pleistocene (Tipple et al., 2010) and the δ 13 C CO2 value of ~ −8‰ corresponding to the time at which the plants were collected (Kohn, 2010). In this paper, we have used the enameldiet enrichment factor of +11.8 ± 0.3‰ calculated for chimpanzees (Pan troglodytes) at Ngogo in Kibale National Park in Uganda for macaques (Malone et al., 2021) and the enamel-diet enrichment factor of +14.1 ± 0.5‰ for ungulates (Cerling and Harris, 1999;Passey et al., 2005) to interpret tooth enamel δ 13 C values and estimate potential consumption of C 3 -or C 4 -plant. These isotopic enrichment factor values for C 3 and C 4 plants are assumed to provide tooth enamel values of < −18.2‰ for closed canopy diets, < −13.7‰ for pure C 3 diets, < −9.7‰ for C 3 -dominated diets, from −9.7 to −0.7‰ for C 3 -C 4 mixed feeders, >−0.7‰ for C 4 -dominated diets and > +1.3‰ for pure C 4 diets in primates and tooth enamel values of < −15.9‰ for closed canopy diets, < −11.4‰ for pure C 3 diets, < −7.4‰ for C 3 -dominated diets, from −7.4 to +1.6‰ for C 3 -C 4 mixed feeders, > +1.6‰ for C 4 -dominated diets and > +3.6‰ for pure C 4 diets in ungulates.
The δ 18 O values in enamel bioapatite is directly related to the drinking water (i.e., springs, rivers, lakes, water holes), plant/food water (i.e., leaves, grass, fruit) and metabolic water by reaction with carbohydrates. Consequently, the precipitation, altitude, humidity/ aridity and local temperatures will have an effect on the δ 18 O of locally available water and food (Longinelli, 1984;Luz and Kolodny, 1985;Bryant and Froelich, 1995;Pederzani and Britton, 2019). Oxygen values can also vary within the same habitat, especially in tropical forest  Crowley, 2014;Fannin and McGraw, 2020). The approximate initial and final mineralization estimated for the extant Macaca nemestrina (Sirianni and Swindler, 1985), in relation to the teeth analyzed in fossil M. cf. sylvanus, is ca. birth to 1 year for M 1 , I 1 ca. 1 to 2 years for I 1 , ca. 1 to 2.5 years for M 2 and ca. 3 to 4 years for M 3 . Related to dental development, the only permanent tooth to begin development prenatally in the non-human primates is the first molar, which may have one to three cusps mineralized at birth and a complete development of the crown for about 1 year (Sirianni and Swindler, 1985). Relatively recent studies in Macaca mulatta suggest that the longest weaning periods are ca. 10 months and the shortest weaning period, between 2 and 5 months (Reitsema et al., 2015). Considering the data on nursing and mineralization timing in extant Macaca, we assume that the I 1 , M 2 and M 3 enamel apatite should reflect a post-weaning signal and the M 1 enamel apatite could therefore reflect a pre-weaning signal for the M. cf. sylvanus from GFT-4.2. This allows us to reconstruct the diet from birth to approximately 4 years of age.
The bivariate graphs were made with the Past 4.02 software (Hammer et al., 2001).

Occlusal microwear analysis (low magnification)
Enamel microwear features were observed via standard light microscopy using a Zeiss Stemi 2000C stereomicroscope at 35× magnification on high-resolution epoxy casts of teeth, following the cleaning, molding, casting, and examination protocol developed by Solounias and Semprebon (2002) and Semprebon et al. (2004). We used the wear facets of the mesiobuccal cusp of M2 (normally the paracone of M2). Before molding, the occlusal surface of each specimen was cleaned using acetone and then 96% ethanol. The surface was molded using high-resolution silicone (Heraeus Kulzer, PROVIL novo Vinylpolysiloxane, Light C.D. 2 regular set), and transparent casts were created using clear epoxy resin (C.T.S. Spain, EPO 150 + K151). A standard 0.16 mm 2 ocular reticle was employed to quantify the number of small pits (round scars-relatively shallow, refract light easily and consequently appear bright and shiny) and large pits (round scarsdeeper, wider, and consequently less refractive), coarse scratches (elongated scars with parallel sides-narrow, and barely etched into the enamel surfaces) and fine scratches (elongated scars with parallel sides-wider and more obviously etched into the enamel surface), gouges (features with irregular edges and much larger and deeper than large pits), puncture pits (craterlike features with regular margins), and cross scratches (oriented perpendicularly to the majority of scratches; Supplementary Figure 1). The scratch width score (SWS) was obtained by giving a score of "0" to teeth with predominantly fine scratches per tooth surface, "1" to those with a mixture of fine and coarse types of textures, and "2" to those with predominantly coarse scratches. We took microphotographs using a Blackfly S digital camera and the Kivy Mic Capture Z software. We used the Helicon Focus 7 software to merge images from various focal planes and to produce a single image with a greater depth of field. We added scale bars using ImageJ. We compared the quantitative results with the reference dataset of the extant Cercopithecidae species. A single observer (IRP) analyzed all specimens to avoid inter-observer errors. We used the Addinsoft™ XLSTAT-3.02 statistical packages to generate the linear discriminant analysis (LDA) and the correspondence analysis (CA).

Buccal microtexture analysis
The original specimens (Supplementary Table 1) were cleaned with acetone and alcohol-soaked cotton swabs, and vinyl impressions were made using President's Jet Regular Body (Coltène-Whaledent).
Positive casts were made with Ferropur (Feroca ® ) polyurethane (Martínez et al., 2022). Buccal areas were obtained using a Sensofar Plu Neox laser scanning confocal microscope, with a 20× (0.45NA) objective, a spatial sampling of 0.83 μm, an optical resolution of 0.31 μm, a vertical resolution of 20 nm, and a z-step interval of 1 μm (Ibáñez et al., 2020;Martínez et al., 2022). Buccal surfaces were placed perpendicular to the objective in a horizontal position and an image of the field of view of the microscope (600 μm × 500 μm) was obtained (see 3D photosimulations in Supplementary Figure 3). From this image, four subsamples of 138 μm × 102 μm (a standard area in DMTA studies) were cropped using the Mountain 7 ® software from Digital Surf. 1 We have tested the S-filter prior to the obtention of the 38 parameters, following (Winkler et al., 2022). The workflow included leveling the surface to remove plane inclination, eliminating the curvature applying a 2nd degree polynomial, removing the anomalous points by applying upper and lower 0.5% (a total of 1% of threshold values), following the threshold operator included in the same software and applying a Spatial filter (Gaussian 3 × 3). The Gaussian filter was applied because is widely used for technical surfaces and was applied after the trimming because it is sensitive to outliers (Schulz et al., 2013). Surface roughness was later quantified according to 38 ISO 25178/ISO 12781 parameters (Martínez et al., 2022;Supplementary Table 5).
To standardize the measured area on the buccal surface, we calculated the median of the parameters derived from the several (one to four) measurements of a single buccal surface (Scott et al., 2006;Calandra et al., 2012;Martínez et al., 2022), and used the median per specimen for analysis. We used log-transformed texture data for analyses because some texture parameters have non-normal distributions (Shapiro-Wilk, p < 0.05). To ensure we did not exclude any parameter (e.g., Ssk has negative values), we log transformed the data as log [1 + parameter] (Supplementary Table 6). We identified the outliers for both original and log-transformed variables according to their Studentized residual values (absolute value > 3) and removed them from the final data set (Schulz et al., 2013;Martínez et al., 2022). To test the null hypothesis that microwear textures do not differ between teeth from different species, we applied analysis of variance (ANOVA) with pairwise testing (Tukey's honest significant difference [HSD]) to each texture parameter for each species. Where homogeneity of variance test (Levene's test) revealed evidence of unequal variances, we used Welch's ANOVA. We selected the parameters that showed significant differences between species and applied a linear discriminant analysis (LDA). Descriptive statistics, tests and graphs were performed using Addinsoft™ XLSTAT-3.02 statistical packages.
Frontiers in Ecology and Evolution 06 frontiersin.org

Occlusal microwear analysis (low magnification)
The M. cf. sylvanus dental microwear pattern is characterized by a high number of pits (NP = 28.5) and a medium number of scratches (NS = 16). The extant cercopithecids show a similar pattern, namely a higher number of pits than scratches (Table 1). Lophocebus albigena has the highest number of pits (NP = 40.6; raw data in Supplementary Table 3 (Table 1).
We performed correspondence analysis to compare four microwear variables (pits, scratches, puncture pits, and gouges) of the extant species and macaque from GFT-4.2 (Figure 1). We plot the results for discriminant functions 1 and 2 (DF1 and DF2, respectively) because they explain a higher percentage of the variance than the other axes (Supplementary Table 4). Three groups are well defined in the graph. The first group is located on the left and comprises T. gelada. The second group, located in the middle, includes M. sphinx. The third group, on the right, includes L. albigena, C. atys, and M. mulatta. T. gelada falls far from the other species because it has the highest number of scratches (NS = 25.9), the smallest number of pits (NP = 10.1), and the absence of puncture pits and gouges. M. sphinx has intermediate microwear values relative to the other two groups. This species has more pits than T. gelada, but also more scratches than the components of the third group. M. sphinx does not have puncture pits. The last group (L. albigena, C. atys, and M. mulatta) is characterized by a high number of pits (NP > 30) and a small number of scratches (NS < 12). The two specimens of the fossil macaque are located near the third group.
We performed LDA including the four microwear variables (pits, scratches, puncture pits, and gouges) to predict group membership of the fossil M. cf. sylvanus specimens from GFT-4.2. The LDA showed significant between-species differences (Wilks' λ = 0.05; F = 20.826; p < 0.0001) and correctly classified 90.48% of the species. The two first discriminant functions explained 96.68% of the total variance. DF1 (85.88%) was positively correlated with scratches and negatively correlated with the presence of puncture pits and pits. DF2 (10.80%) was Frontiers in Ecology and Evolution 07 frontiersin.org Correspondence analysis based on four occlusal microwear variables (number of puncture pits, number of gouges, number of scratches, and number of pits). Ellipses show 68% confidence limits for means (one standard deviation). Eigenvalues, percentages of variance and cumulative percentages of each axis in Supplementary negatively correlated with the presence of gouges. GFT140 was classified as C. atys (99% of predicted classification probability) and GFT138 was classified as L. albigena (83% of predicted probability).

Buccal microtexture analysis
A total of 15 buccal enamel microtexture parameters out of the 38 (39.47%) showed significant differences (ANOVA and Tukey's HSD post hoc procedure p < 0.05 and Welch correction, see Supplementary Table 7) among the Cercopithecoidea species. The main significant differences are between the folivore-frugivore M. mulatta and the graminivorous T. gelada and the hard fruit eater C. atys and the folivore-frugivore M. sphinx. T. gelada mainly differentiate from M. mulatta in volume parameters, functional and flatness (Vv, Vmc, Vvc, Sxp, Smc, FLTp, FLTv, and FLTq p < 0.05) while the hard fruit eater differentiates from M. sphinx in spatial, feature and flatness parameters (Sal,Shv,FLTp,FLTt,FLTv,and FLTq;Supplementary Figure 4). No significant differences are found between hard fruit eaters taxa, the arboreal L. albigena and terrestrial C. atys and, only one difference was found between the hard fruit eaters and the graminivorous T. gelada (the feature parameter Spc between L. albigena and T. gelada p < 0.05).
The LDA includes the 15 parameters to classify the extinct GFT-4.2 samples. DF1 and DF2 account for the 82.51% of the total original variance (Wilks' λ = 0.014; F = 1.371; p-value = 0.143). DF1 (57.15%) is positively correlated with flatness and volume parameters and negatively with Shv (mean hill volume) DF2 (25.36%) was positively correlated with all the parameters except for Spc (Figure 2; see Supplementary Table 8  The mean heigh, volume and flatness parameters of the three specimens from GFT-4.2 differ significantly (p < 0.05) from both C. atys and T. gelada, and were similar to folivores and seed eaters. Two GFT-4.2 specimens, therefore, resemble to primates that exploit seasonally leaves, fruits and seeds in major proportion, on the ground of mixed forest. GFT137, however, is classified as T. gelada a graminivorous primate that exploit also forbs in open ecosystems, due to similarities on buccal microtexture feature parameters.

Discussion
Our multiproxy analysis suggests a different short-term (occlusal) and long-term (buccal) feeding habits of M. cf. sylvanus from GFT-4.2, but with an overall dominance of abrasive, tough and brittle food items in an environment dominated by C 3 plants. The most commonly consumed foods would be fruit, seeds and leaves but also some grasses (Table 2). This dietary variety reflects the ability and flexibility of fossil macaques to consume different foods. As well as the richness of resources, and probably microhabitats, present in the GFT-4.2 area during the Plio-Pleistocene.
Based on short-term resolution, the occlusal microwear analysis classifies the M. cf. sylvanus specimens (GFT138 and GFT140) as L. albigena and C. atys, which have a diet rich in hard fruit and seeds (Poulsen et al., 2001;Shah, 2003;Lambert et al., 2004;Daegling et al., 2011;McGraw et al., 2011McGraw et al., , 2012McGraw et al., , 2014. The extant L. albigena and C. atys species and the fossil M. cf. sylvanus have puncture pits and gouges, which is a distinctive microwear feature of hard-object feeders Semprebon et al., 2004;Williams and Geissler, 2014). Both Cercocebus and Lophocebus have thicker enamel compared to other primate groups, which allows both genera to rely on hard foods throughout the year (McGraw et al., 2012). Macaca cf. sylvanus occlusal surfaces differs significantly from the extant T. gelada, which is characterized by an almost exclusively graminivorous diet with seasonal incorporation of forbs (Fashing et al., 2014;Jarvey et al., 2018). The occlusal microwear pattern of T. gelada is dominated by a high number of scratches and the absence of puncture pits and gouges. This surface is characteristic of animals with grass-dominated feeding habits Semprebon et al., 2004). Even if it was not the dominant food item in their diet, we cannot exclude that M. cf. sylvanus consumed Plot of the two first discriminant functions (DF1-2) of buccal-dental microtextural parameters. Ellipses show 68% confidence limits for means (one standard deviation). On the top left are the loadings of the ISO texture parameters onto DFs. Raw data in Supplementary Table 8 (illustrations by Diego Rodríguez-Robredo).
Frontiers in Ecology and Evolution 09 frontiersin.org some grasses due to the large number of scratches found in specimen GFT138. Occlusal microwear and buccal microtexture analysis group fossil M. cf. sylvanus with different extant species, an outcome that is explained by differences in the temporal resolution and turnover of the analyzed tooth surfaces. Buccal microtexture analysis provides medium longterm, compared to occlusal microwear, from weeks to months (Romero et al., 2012) temporal resolution classifying two of the three M. cf. sylvanus specimens (GFT 138, and GFT 140) from GFT-4.2 as M. sphinx and M. mulatta. These species primarily consume fruits, leaves, and seeds in closed and riparian forests during the wet season, and seeds foraged from the ground year-round. The sample included in the present analysis comes from the population that inhabit Pakistan where M. mulatta rely primarily on leaves and fruits depending on their availability (Goldstein and Richard, 1989;Sengupta and Radhakrishna, 2016). Both taxa, therefore, consume abrasive items from the ground in mixed forest. Mandrillus sphinx and M. mulatta diet also includes herbs such as piths and blades and occasionally woody tissue during seasons of heightened fruit scarcity (Sengupta and Radhakrishna, 2016;Hongo et al., 2017). Based on observations of M. sphinx from Gabon, more seeds are consumed during the dry season than in the peak fruit season. Moreover, finely crushed seeds in fecal samples indicate that M. sphinx uses the occlusal surface to fracture seeds (Hongo et al., 2017), and this mechanical action correlates with gouges on these surfaces (Figure 3;  Supplementary Figure 2). M. sphinx consumes large nuts, like S. gabonensis, during seasons of heightened fruit scarcity (Hongo et al., 2017). Hard items are consumed habitually by L. albigena and C. atys year-round (McGraw et al., 2012). The microtexture buccal pattern of M. sphinx, characterized by higher relief, with peaks and dales and coarse features compared with the rest of analyzed species, may correlate with the presence of more abrasive fractured seeds or phytoliths from the leaves that are processed in the oral cavity while being fractured along occlusal surfaces. The similarity in the microtexture pattern of the two GFT-4.2 specimens with M. sphinx and M. mulatta suggests that M. cf. sylvanus regularly consumed abrasives food items. The LDA of the buccal microtexture plot two specimens (GFT138 and GFT140) close to each other and to the hard fruit eaters (Figures 1,2). All these species rely on hard-shelled fruits and seeds with seasonal fluctuations. Specially, M. mulatta diet resembles M. sphinx (Sengupta and Radhakrishna, 2016) with more soft fruit ingestion in the former. Both M. mulatta and M. sphinx incorporate hard and brittle dietary resources but, this is dependent on fruit availability during the year (Sengupta and Radhakrishna, 2016). When fruit is scarce, M. mulatta relies heavily on seeds. However, the buccal microtexture pattern of the GFT137 is also similar to T. gelada (Figure 2). Theropithecus gelada is a graminivorous primate that inhabit the open highlands of Ethiopia and consume grasses and forbs (Jarvey et al., 2018;Souron, 2018).
Three out of four extant species in which the microwear classifies the macaque of GFT-4.2 are species whose diet is related to the presence of trees. In this sense, oxygen values could also reinforce this link. The four M. cf. sylvanus specimens (GFT137, GFT138, GFT139, and GFT140) from GFT-4.2 show high δ 18 O values compared with the contemporary ungulates recovered from the same site. The same pattern occurs with other fossil cercopithecines (T. atlanticus and T. oswaldi) from North Africa (Figure 4). Extant M. sylvanus are highly adaptable and have an eclectic diet which adds to the variables that can affect their oxygen isotope values. A significant amount of their water ingestion comes from leaf, fruit and vegetable consumption, especially during the rainy season. Meanwhile, during the dry season, they drink more water from open water sources (Fooden, 2007). In the case of M. cf. sylvanus TABLE 2 Summary of the results provided by the three methods used to reconstruct the diet and habitat of Macaca cf. sylvanus from GFT-4.2.

Methods
Diet and habitat reconstruction Timescales Buccal microtexture Folivore-frugivore with hard item consumption and some grass intake (Similar to extant Mandrillus sphinx, Macaca mulatta, and

Theropithecus gelada).
Toward the end of the animal's life (Long-term/weeks to months) Occlusal microwear Hard-fruit eater with leaf and seed intake (Similar to extant

Lophocebus albigena and Cercocebus atys).
Toward the end of the animal's life (short-term/last days) The timescale represented by the enamel isotopic signature has been defined following the work of dental development in the species Macaca nemestina (Sirianni and Swindler, 1985) and the weaning period defined for M. mulatta (Reitsema et al., 2015). Frontiers in Ecology and Evolution 10 frontiersin.org from GFT-4.2, the differences in oxygen values in relation to contemporary ungulates could be associated with differential access to and consumption of certain food resources by the primates, such as fresh fruits and leaves from trees. Sun-exposed leaves, especially in arid conditions, are particularly sensitive to evaporative fractionation that could lead to 18 O enrichment (Gonfiantini et al., 1965;Barbour, 2007;Crowley et al., 2015;Carter and Bradbury, 2016;Fannin and McGraw, 2020). Other authors argue that this enrichment by evapotranspiration also occurs in other parts of the plant like fruits (Yakir, 1997;McCarroll and Loader, 2004). Paleoenvironmental data of GFT-4.2 based on the study of paleoherpetological and small mammal records, suggest dry conditions in an open woodland habitat with a permanent water body nearby Piñero et al., 2019;Alba et al., 2021).
Considering that aridity affects all GFT-4.2 taxa, we suggest that the oxygen enrichment of M. cf. sylvanus may be due to the ingestion of some different foods. An example in extant primates is the one exposed by Moritz et al. (2012) that measured the δ 18 O hair from two baboon species (Papio anubis and Papio hamadryas) which inhabit the arid hybrid zone of Awash National Park (Ethiopia). The diets of the two species in this area are similar and both obtain their meteoric water from the same source (Awash River). However, the oxygen values of P. hamadryas are relatively enriched due to the ingestion of a greater proportion of 18 O-enriched leaf water in the arid thornbush. Krigbaum et al. (2013) show how δ 18 O ap values of different cercopithecids covary depending on foraging position in the closed tropical canopy forest (Taï forest, Côte d'Ivoire). Higher humidity on the forest floor is associated with low δ 18 O ap values, with enrichment of 18 O observed along a vertical gradient. The C. atys exhibits low δ 18 O ap values relative to guenons and colobines, which reflects the C. atys terrestrial foraging habits. Although we do not have a tropical canopy forest in GFT-4.2, these differences in δ 18 O values between feeding on the ground and in the upper parts of the tree could explain the 18 O enrichment of macaques compared to ungulates.
The 18 O enrichment points to a diet based on tree foraging in a forest ecosystem. These findings are consistent with the microtexture and low magnification signals and support the carbon isotope results of M. cf. sylvanus from GFT-4.2. The δ 13 C values indicate consumption of C 3 plants (i.e., trees, shrubs, fruit/seeds, and C 3 -grass adapted to a mild growth season). M. cf. sylvanus, Hipparion sp., Gazella sp. and Tragelaphini from GFT-4.2 have been compared with the fossil T. atlanticus from Ahl al Oughlam, Morocco (ca. 2.5 Ma) and T. oswaldi from Tighenif, Algeria (ca. 1 Ma) and other ungulates from the same sites (Fannin et al., 2021;Figure 5;Supplementary Table 2). The M. cf. sylvanus (n = 4) mean is −12.3 ± 0.3‰, the T. atlanticus (n = 10) mean is −12.4 ± 0.8‰, and the T. oswaldi (n = 11) mean is −11.3 ± 1.9‰ ( Figure 5). These three fossil species show similar median carbon values and no consumption of C 4 plants. Except for one individual of T. oswaldi from Tighenif with a δ 13 C value of −5.9‰ that includes some C 4 plant consumption.
Frontiers in Ecology and Evolution 12 frontiersin.org indicate a preference toward food sources associated with closed woodland environments (Adams et al., 2013). This interpretation is consistent with the carbon values of all of them, including the fossil macaque. Feeding habits of P. robinsoni and P. jonesi from Swartkrans have been reconstructed as a frugivore-folivore on the basis of cusp form (Benefit, 2000). We also propose this using the occlusal and buccal microwear for M. cf. sylvanus from GFT-4.2, although the exploitation of grasses and forbs in a more open habitats cannot rule out due to the buccal microtexture pattern of GFT137 and the high number of scratches on the occlusal surface of GFT138. Carbon isotopic signal of the Macaca from GFT-4.2 would indicate that these consumed grasses would be C 3 . Consumption of C 3 grasses has also been identified in the diet of T. oswaldi from the Plio-Pleistocene site of Ahl al Oughlam (Morocco), as well as in extant T. gelada at Guassa Plateau (Ethiopia; Fannin et al., 2021). The microwear analysis shows that the diet of the fossil macaque is close, among others, to that of extant M. sphinx, L. albigena, and C. atys, species that live in dense tropical forests. Isotopic studies on extant M. sphinx and C. atys have shown carbon values strongly affected by canopy effect. According to Oelze et al. (2020), the δ 13 C values of M. sphinx hair sections from Gabon are within the range (−25.1‰ ± 0.3‰) of those obtained from the hair of chimpanzees and gorillas in Loango National Park. Like these species, mandrills live in dense forest habitats. Another study by Krigbaum et al. (2013) on bone apatite from seven sympatric cercopithecine species, including C. atys, in the Taï forest (Côte d'Ivoire) obtained carbon values range between −19.9‰ and −17.9‰ for all species. In contrast, the carbon values of M. cf. sylvanus from GFT-4.2 (−12.3 ± 0.3‰) do not correspond to those of a species foraging in a closed canopy forest. Considering the enamel-diet enrichment factor for primates (+11.8 ± 0.3‰; Malone et al., 2021) and the cutoff of −30‰ (δ 13 C CO2 corrected by ~ +1.5) for closed-canopy forests (Kohn, 2010), the δ 13 C values of tooth enamel for a diet in a closed-canopy forest should be lower than ~ −18‰ for fossil Macaca. The carbon values of the taxa analyzed are compatible with an open woodland. This is in line with the paleoherpetological and small mammal work in which taxa characteristic of open woodland have been defined Piñero et al., 2019).
The use of these three methods has allowed us to produce a paleodietary reconstruction at different times in the life of the individuals, although we cannot discard that the analyzed molars belong to the same specimen. The isotopic signal, corresponding to approximately the first 4 years of life, is in agreement with the data of the two dental microwear proxies in the diet-tree relationship for M. cf. sylvanus, and the absence of C 4 typical of open grassland. Buccal and occlusal microwear turnover is a dynamic process; however, the buccal signal is more stable in the long term and may indicate the diet consumed year-round (Romero et al., 2012). Whereas the low magnification occlusal microwear results may indicate the consumption of hard fruit and seeds similar to mangabeys the last days or weeks prior to death , perhaps representing the diet of a single season. The buccal microtexture analysis suggest that the abrasiveness of their diet was similar to primates (M. sphinx and M. mulatta) that consumes hard, tough, and brittle food items year round with the incorporation of leaves, and the consumption of grasses as T. gelada does. The complementarity of these proxies helps to differentiate between seasonal, short-term durophagous species from those that exploit tough and brittle items actively in a long term, broadening the dietary range of consumed resources throughout the year.

Conclusion
Our research represents the first report related to the feeding ecology of fossil Macaca during the Plio-Pleistocene in Africa. Our results highlight a tight correlation between the fossil M. cf. sylvanus specimens from GFT-4.2 and extant primates with a predominant abrasive, tough and brittle food items. Evidence from the multiproxy analysis suggests that the macaques foraged for food in a wooded habitat, but also foraging in more open environments. The versatility of macaques to occupy wide habitat ranges makes them a useful candidate for analogies with hominins and for understanding how they adapted to different environments. Evidence of early hominins in North Africa is dated close to the Plio-Pleistocene boundary (ca. 2.4 Ma) at the Ain Boucherit site in Algeria. The ecological context of this first human population is a key question to understand the dispersions of humans and other mammals in these territories. Considering the proximity of GFT to Ain Boucherit, understanding the feeding ecology of M. cf. sylvanus from GFT may provide clues about some of the ecological resources that these early humans may have had. This initial investigation is only one step toward understanding the feeding habits of M. cf. sylvanus. More data are needed from additional individuals to be able to confidently interpret foraging ecology of this species.

Data availability statement
The original contributions presented in the study are included in the article/Supplementary material, further inquiries can be directed to the corresponding authors. The surface data are available from the corresponding author upon reasonable request.  AR-H, and CT research was funded by CERCA Programme Generalitat de Catalunya. IR-P is beneficiary of predoctoral fellowship (2020-FI-B-00731) funded by AGAUR and the Fons Social Europeu (FSE). AA and is beneficiary of a fellowship from the Erasmus Mundus Program to do the Master in Quaternary and Prehistory at the Universitat Rovira i Virgili (Tarragona, Spain). CT was supported by the Spanish Ministry of Science and Innovation through the "Ramón y Cajal" program (RYC2020-029404-I). The Institut Català de Paleoecologia Humana i Evolució Social (IPHES-CERCA) has received financial support from the Spanish Ministry of Science and Innovation through the "María de Maeztu" program for Units of Excellence (CEX2019-000945-M), including the postdoctoral fellowships of AR-H.

Acknowledgments
We thank the Direction of Cultural Heritage (Ministry of Culture and Communication, Morocco), the Jerada Government and Local Authorities of Aïn Béni Mathar and Guefaït for local permits to develop geological, archeological and paleontological fieldwork in the region. Lastly, we want to thank the local population for their support, collaboration and help in the discovery of some of the fossils and certain archeo-paleontological sites and above all for their hospitality and welcome. We also acknowledge the two reviewers, the editor Brooke E. Crowley and Robert Patalano for the constructive comments raised, which considerably helped to improve this research and Diego Rodríguez-Robredo for the illustrations of extant primates. We are grateful to all the curators and technical personnel of the different institutions where the extant specimens were molded.

Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher's note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.