Diet of the Insular Lizard, Podarcis lilfordi (Günther, 1874): Complementary Morphological and Molecular Approaches

Simple Summary Traditionally, the diets of lizards and other small vertebrates have been studied using invasive techniques, such as the analysis of gastric contents obtained by dissection or stomach-flushing. Nowadays, the morphological analysis of the remains contained in feces is commonly used. However, these techniques require a great deal of experience to identify prey remains, and some prey items may be left undetected in the analyses. Recently, the use of molecular techniques has made it possible to identify prey from feces, thus significantly expanding the diversity of identified prey. Despite this, the molecular analysis of diets also has limitations, since it does not allow an adequate quantification of the contribution of each prey type, and, on occasion, some prey may not leave traces of DNA detectable in the feces. For this reason, here we propose a study of the diet in which the two methods of fecal analysis are complementarily used: morphological identification and identification by means of DNA analysis. Abstract The diets of insular lizards are extremely varied, depending on the different environmental characteristics of each island population. This is particularly evident in the case of the populations of small coastal islets of the Balearic Islands, where the Balearic lizard, Podarcis lilfordi, is found. The study of trophic ecology carried out by means of traditional tools, such as morphological analysis of feces, has made it possible to detect numerous prey and nutritional elements. However, these methods are clearly insufficient, as some rare groups are not detected. It is also difficult to identify remains of marine subsidies or of foods contributed to these small islands by other predators, such as seabirds. The current study demonstrates the advantages of combining morphological diet analysis with the molecular study of individual feces samples obtained from the same populations. We obtained a greater diversity of prey groups using the combined methodologies, with each method identifying prey items that were not detected using the other method. Particularly, the study of diets at the molecular level identified plant species consumed by lizards that were, occasionally, not identified in morphological analyses. Conversely, the traditional morphological study of an equivalent number of fecal samples allowed for the identification of several prey groups that had not been detected in the molecular study. From this viewpoint, the advantages and disadvantages of each methodology are discussed.


Introduction
Lizard diets can be complex, and feeding habits are expected to vary among species depending on evolutionary history, microhabitat characteristics, and prey availability [1]. Most lizard species are generalist predators [2]. Traditionally, they were considered insectivores, although collectively they consume a much wider variety of food types, including stored in absolute ethanol vials. Lizards were immediately released at the point of capture. The samples were preserved at 4°C in the field; upon arrival to the laboratory, they were stored at −20 °C until DNA extraction [12]. Additionally, a total of 138 different individual feces samples, from the same season and localities as the molecular samples, were included in the morphological analysis. For these samples, we collected feces from individual lizards, as well as from the ground.  Additionally, a total of 138 different individual feces samples, from the same season and localities as the molecular samples, were included in the morphological analysis. For these samples, we collected feces from individual lizards, as well as from the ground. to the Roy J. Carver Biotechnology Center (University of Illinois, Ill., USA) for amplification in a microfluidic high-throughput multiplexed PCR platform (Fluidigm). For animal prey detection, we used mitochondrial cytochrome oxidase I (COI) for the following primer pairs: mlCOlintF/jgHCO2198 (5 -GGWACWGGWTGAACWGTWTAYCCYCC-3 /5 -TANACYTCNGGRTGNCCRAARAAYCA-3 [20,21] and specific arthropod ArtF11/ ArtR17 (5 -GGNKYNGGNACWGGATGAACWGTNTAYCCNCC-3 /5 -GGRTCAAAAAA TGAWGTATTHARATTTCGRTCWGTTA-3 [22]. To study the ingested plants, the following rbcL and psbA-trnH chloroplast markers were used: psbA3_f/trnHf_05 (5 -GTTAT GCATGAACGTAATGCTC-3 /5 -CGCGCATGGTGGATTCACAATCC-3 ; [23,24] and rb-cLa_F/rbcLa_R (5 -ATGTCACCACAAACAGAGACTAAAGC-3 /5 -GTAAAATCAAGTC CACCRCG-3 [25,26]. CS1 and CS2 Fluidigm universal tags and barcode labels specific to each sample and Illumina adapters i5 and i7 were used. The resulting amplicons were validated on a Fragment Analyzer (Agilent) using the HS NGS kit (DNF-474-33). Sequencing was conducted on an Illumina MiSeq v2 platform yielding 2 × 250 paired-end reads. Regarding animal prey, a blocking primer-targeting region within COI was used to minimize host amplification. For the design, we used COI sequences of the most commonly consumed prey reported from previous studies of Balearic lizards [9] available at Barcode of Life Data System (BOLD; http://www.boldsystems.org, accessed on 1 September 2022) and COI fragments of P. lilfordi available through the GenBank accession codes in [18].

Sequence Analyses and Taxonomic Assignment
Micca version 1.7.2 [27] was used for merging, trimming, filtering, and OTU picking (see [12] for more details). For the taxonomic assignment (Table 1), the resulting DNA sequences (OTUs) were joined into 98% similarity clusters with Usearch version 10.0.240_i86osx32 [28]. Consensus sequences (centroid) were defined for each cluster and were used to download, from public DNA databases (GenBank and rrnDB), the 1000 most similar sequences.
Individual centroids were analyzed using all the sequences forming the cluster and the 1000 most similar sequences. DNA matrices were aligned with MAFFT [29] and used for phylogenetic inference employing IQTREE [30]. Trees were explored in FigTree [31] to establish the systematic position of the diet sequence to the highest taxonomic rank according to bootstrap support values ≥ 70% [32]. Information from multiple markers targeting the same group of organisms (i.e., animals or plants) was merged following the multi-marker metabarcoding approach described in [33].
Finally, taxonomic assignments and current nomenclature from metabarcoding (Table 1) were checked with available information provided in the taxonomic literature of plants and invertebrates recorded in the Balearic Islands (List A1 (prey references) and List A2 (plant references) from Appendix A).

Morphological Study of the Diet
We analyzed fecal pellets under a binocular dissecting microscope. Diet reconstruction was based on a meticulous pellet analysis that is highly comparable to diet reconstructions based on gastric contents, with soft-bodied prey and, particularly, insect larvae being equally represented in fecal pellets and gut contents [10]. Each scat was spread in a thin layer of less than 0.5 mm over the entire surface of a Petri dish with some drops of 70% ethanol. Plant families were arranged according to Cronquist [34] and Bremer et al. [35]. Prey remains were visually identified up to their order or, in some cases, family level (Table 2). Prey number for each fecal pellet was conservatively estimated by counting only easily identifiable remains. We detected the consumption of plant remains in morphological analysis, but without any taxonomic assignment (see [9] for more details).

Diet Comparison
Because fecal pellet sample sizes were different between morphological and molecular analyses, we tested for the contribution of each prey type to the diet between molecular and Animals 2023, 13, 507 5 of 23 morphological samples only, using a subset of 10 feces from larger morphological samples (Table 2). This subset was randomly extracted from each population sample. Because we could not estimate the quantitative contribution of each prey type from molecular samples, we compared only their presence [9] in the feces.
For the five populations, the presence of different prey types found with morphological and molecular methods (Table 2) of analysis were compared with a permutational multivariate analysis of variance (permutational MANOVA, [36]), using the 'adonis' function from 'vegan' R package [37]. The multivariate homogeneity of group dispersions (variances) was tested with the function 'betadisper', a multivariate analogue of Levene's test for homogeneity of variances.   We estimated and compared diet diversities using the approach proposed by Pallmann et al. [38]. We converted "raw" Simpson or Shannon indices into "true" diversities, which belong to the same mathematical family. Thus, different measures of area were regarded as special cases of Hill's general definition of diversity measures [39]. Two-tailed tests for integral Hill numbers were performed to compare diets obtained from the different analytical methods. This selection includes the transformed versions of the two following indices: the Shannon entropy index, H sh (q→1), and the Simpson concentration index, H si (q = 2, [40]). We performed 5000 bootstrap replications to obtain reliable p-values [41]. Methods described here were implemented in the R package "simboot" [42] and are fully described in [38]; more details of dietary analysis can be found in [9].
The frequency of each prey type present in the fecal samples was compared with the Fisher exact tests. All calculations were performed in R version 4.0.3 [43].

Morphological and Molecular Diet Compositions
In molecular analysis, a total of 57 prey taxa were identified corresponding to 17/5 orders and 27/6 families (Animals/Plants, respectively) ( Table 1). The most common orders iden-Animals 2023, 13, 507 7 of 23 tified for animal components were Isopoda, Hymenoptera, Lepidoptera, and Coleoptera ( Figure 2). Within order Hymenoptera, Formicidae was the only represented family, while Halophilosciidae was the most abundant of the Isopoda orders ( Figure 3). In morphological analysis, we identified 9 different prey taxa in the random sample of 10 fecal pellets (Table 2). With regards to plant components, Caryophyllales ( Figure A1) was the most abundant Order, represented by Amaranthaceae and Plumbaginaceae families ( Figure 3).

Comparison of Morphological and Molecular Analytical Methods
Comparisons were made, taking into account only ten randomly-selected morphological molecular and molecular samples. Results are represented in Figure 2.
We found overall significant differences in diets from morphological and molecular analyses (adonis, F 1,79 = 6.7728, p = 0.0009, with homogeneous variances, F 1,79 = 0.022, p = 0.88). However, we obtained variable results from the comparison of subsamples from each population. At the Ses Bledes islet, we did not find significant differences in the presence of different taxa between the two samples of ten fecal pellets (Fisher Exact Test, p = 0.3038), because the diet comprised mainly ants (Tables 1 and 2). Regardless, the detection of different prey items, using the two methodologies, was extremely different, with four groups of prey not detected by molecular analysis (Gastropoda, Pseudoscorpionida, Hymenoptera other than Formicidae, and carrion from dead birds) and three groups undetected in the morphological random subsamples (Araneae, Diptera, and Lepidoptera). A similar situation was observed in Caragol (Fisher Exact Test, p = 0.6233), where the diet was dominated by Isopoda, Coleoptera, and Formicidae (Table 1) and in Na Guardia, where a very diverse diet was detected, with a small presence of each prey type in both samples (Fisher Exact Test, p = 0.2841). In the case of the Moltona Islet, molecular and morphological results were very different, because of a higher detection of groups using the molecular approach (Fisher Exact Test, p = 0.01125). This was also the case for Na Pelada Islet samples (Fisher Exact Test, p = 6.86 × 10 −5 ).

Comparison of Morphological and Molecular Analytical Methods
Comparisons were made, taking into account only ten randomly-selected mor logical molecular and molecular samples. Results are represented in Figure 2.  Trophic diversity was significantly higher in molecular samples (Shannon index for morphological samples, H sh = 2.239 ± 0.011, and for molecular samples, H sh = 2.451 ± 0.006, p < 0.001; Simpson index for morphological samples, H si = 0.837 ± 0.0005, and for molecular samples, H si = 0.878 ± 0.0002, p < 0.001).

Discussion
Previous literature suggests that the use of metabarcoding allows for greater precision in prey identification, as well as the taxonomic identification of plant elements in the organism's diet [12]. Moreover, a molecular approach may lead to an increase in the number of identifiable prey items [12,44]. In our study, numerous prey items were identified at the Species level, while others were only identified at the Order or Family level using morphological methods of diet analysis (see also [12]). This was particularly true for some groups, such as Diptera, Hemiptera, Lepidoptera, and Coleoptera (Tables 1 and 2). Moreover, molecular analysis allowed the discovery of groups not detected using morphological analysis, such as Psocoptera, Embioptera, or Neuroptera. We showed that, overall, diversities determined by molecular techniques were significantly higher than those determined using similar sample sizes of fecal pellets.
Nevertheless, the study of diet using molecular approaches also holds some limitations and disadvantages, particularly with respect to traditional morphological methods. For example, the use of metabarcoding prevents the identification of cannibalism [44], which occurs frequently in small Mediterranean islets [9]. Additionally, current molecular techniques do not yet allow us to establish the quantitative contribution to the diet of each prey type [44]. They only provide us with a snapshot of the diet of a species-we cannot estimate its diversity, biomass contribution, prey size distribution, and other dietary descriptors.
Unlike other studies, we did not assemble a reference collection of plant species and potential prey in the populations under study; we also note that these identifications are not always reliable. Identified taxa sometimes belong to Species or Genera that are not present in the western Mediterranean or the Balearic Islands. Consequently, a post hoc analysis of these preliminary identifications was necessary, using taxonomic literature. In addition, an increase in the number of prey taxa identified through metabarcoding did not always occur. In most populations, we identified some additional prey using molecular methods, while other prey, even at the order level, were only identified through morphological observation of fecal pellets, and did not appear in the molecular samples. Hence, it is clear that the ideal scenario is to employ both methodologies as complementary approaches to describe the diet of a species [44].
Some molecular findings were difficult to interpret. In the Na Moltona Islet, the presence of DNA from the genus Mus, was detected in lizard scats. The presence of the House mouse (or the Algerian mouse, Mus spretus) on this coastal islet with a long history of human presence [45] cannot be ruled out. In fact, these species have been reported from islands with larger surface areas in the Balearics, but could be absent from most of the coastal islets. If this is the situation, we can posit that DNA from this small mammal was obtained from pellets of such scavengers as the Yellow-footed gull, Larus michahellis (see similar results in [12]).
Other results of the molecular analysis deserve our attention, such as the frequent consumption of several species of crustaceans. In many cases, these are terrestrial isopods, well-known from the diet of P. lilfordi, and present in morphologically analyzed samples [9]. Such is the case with the genus Porcellio. This is also the case with Ligia italica, a very common marine isopod in the diet of this species of lacertid lizards (see, for example, [9]). The significance of these findings, apart from the specific identification, is the confirmation of lizard foraging in shore areas and the frequent consumption of marine subsidies [12].
In the Caragol islet sample, we detected the presence of DNA of true marine crustaceans, particularly Pachygrapsus marmoratus, an extremely common decapod from the shore waters of the Balearic Islands [46], which had already been detected as carrion food of P. lilfordi from Aire Island [9]. This foraging behavior was confirmed by the presence of Cymodoceaceae DNA in samples from the same islet (Table 1). The only species of this seagrass present on the coast of Mallorca is Cymodocea nodosa [47]. In fact, some individual lizards seem prone to forage in shore areas. In Caragol Islet, from a single molecular fecal sample, we detected the presence of P. marmoratus, L. italica, and the seagrass C. nodosa (Table 1).
Regarding the consumption of some plant species, the molecular results must also be interpreted with caution and contrasted with morphological and ecological data. The supposed consumption of the genus Pinus has been detected in three populations where there are no pines, although these populations are close to the coast of Mallorca Island, an area with abundant pine forests. Thus, access to the fruit seems unlikely and the needles do not seem appropriate as a nutritive element [9]. This result poses a problem of a more general nature, i.e., the possibility that the detection of DNA from some plant species may be due to the adherence of pollen grains to the outer surface of some fecal pellets. In our samples, this possibility can be ruled out since the droppings were obtained directly from captured individuals. However, pollen grains of plant species not directly consumed by lizards may be present on the outer surfaces of other prey items [12]. It is important to note that the metabarcoding of plants cannot distinguish between the different parts of a plant species, such as its leaves, pollen, or fruit.
On Guardia and Pelada islets, DNA from the genus Prunus, which includes several species of fruit trees, was detected in feces. There are no species of this genus on either of the two islets. We suspect that its presence in our samples of lizard feces was due to the provision of fruit by visitors to the islets.

Conclusions
We can conclude that the application of molecular techniques in the study of lizard diets is extremely useful in combination with the traditional morphological analysis of feces. Using molecular evidence, we confirmed some direct observations of food consumption and found new taxa not detected using morphological analysis. Moreover, we were able to establish the specific identification of carrion, marine subsidies, and several common preys, such as ants. This study has the potential to be instructive by drawing the attention to new food items [48].

Conflicts of Interest:
The authors declare no conflict of interest.
Appendix A Figure A1. Frequencies of occurrence of each order of plants in the diet composition of Podarcis lilfordi through molecular analysis.
Lists A1 and A2. Additional list of references employed to check the validity of molecular identifications from metabarcoding, according to the taxonomic literature of prey and plants of Balearic Islands.