DOI: 10.5061/dryad.z8w9ghxmm
title: With or without you: Gut microbiota does not predict aggregation behavior in European earwig females

description: The reasons why some individuals are solitary and others gregarious are
 the subject of ongoing debate as we seek to understand the emergence of
 sociality. Recent studies suggest that the expression of aggregation
 behaviors may be linked to the gut microbiota of the host. Here, we tested
 this hypothesis in females of the European earwig. This insect is ideal
 for addressing this question, as adults both naturally vary in the degree
 to which they live in groups and show inter-individual variation in their
 gut microbial communities. We video-tracked 320 field-sampled females to
 quantify their natural variation in aggregation and then tested whether
 the most and least gregarious females had different gut microbiota. We
 also compared the general activity, boldness, body size, and body
 condition of these females and examined the association between each of
 these traits and the gut microbiota. Contrary to our predictions, we found
 no difference in the gut microbiota between the most and least gregarious
 females. There was also no difference in activity, boldness, and body
 condition between these two types of females. Independent of aggregation,
 gut microbiota was overall associated with female body condition, but not
 with any of our other measurements. Overall, these results demonstrate
 that a host's gut microbiota is not necessarily a major driver or a
 consequence of aggregation behavior in species with inter-individual
 variation in group living and call for future studies to investigate the
 determinants and role of gut microbiota in earwigs.
description: <em><strong>Animal sampling and laboratory
 rearing</strong></em> We field sampled
 1000 <em>F. auricularia</em> females and males in a peach
 orchard near Valence, France (Lat 44.9772790, Long 4.9286990) at the end
 of June 2022. All these individuals were sampled as mated adults and
 belonged to <em>F. auricularia</em> species “A” [1, 2]. Upon
 collection, the 1000 individuals were first mixed in a single plastic
 container and then randomly distributed among four other plastic
 containers of 100 females and 100 males each (called “test containers”)
 and two additional plastic containers of 100 females each (called
 “attraction containers”). The test container held females for which we
 then measured the aggregation level, general activity, boldness, fresh
 weight, and gut microbial diversity (see below). The attraction containers
 held females that served as stimuli for measuring aggregation levels. All
 plastic containers were grounded with moistened sand and cardboard
 shelters. Throughout the experiment, the animals were fed with carrot
 chunks <em>ad libitum</em>, which were changed twice a week.
 Containers were maintained under a controlled 12:12 light-dark cycle at
 20°C and 18°C, respectively [3]. All measurements proceeded three days
 after laboratory conditions had been set up.
 <strong>Aggregation measurements</strong>
 To investigate the association between the level of aggregation
 of a female and her gut microbial diversity, we first measured the
 aggregation level of 320 females from the “test containers” (Fig S1). To
 this end, we used 3D-printed arenas consisting of four linearly aligned
 circular chambers (diameter 4 cm) (Fig S1; Van Meyel et al., 2022). Three
 of the chambers were connected by 0.5 cm wide corridors allowing earwigs
 to move between chambers. The fourth chamber (called the “attraction
 chamber”) had a reduced corridor of 0.15 cm to prevent animals from
 reaching the other chamber while allowing the circulation of odors and
 antennal contacts between individuals on both sides [5].
 The measurement of each aggregation score started by placing a
 randomly selected female from the two test containers in the circular
 chambers at one end of the system, and three females taken randomly from
 the two attraction containers in the attraction chamber at the other end
 of the system (we alternated the orientation of the system between
 trials). All these females had access to the food source until they were
 used in the 3D-printed arenas. To avoid injuries during animal handling,
 all these individuals were previously anaesthetized with
 CO<sub>2</sub> for 30 s. One hour later (to allow recovery
 after anesthesia), we started to record whether the focal female was in
 the chamber next to the group of females (yes or no) and then repeated
 this measurement by taking pictures every hour for 72 hours using infrared
 cameras and the software pylon Viewer v5.1.0 (Basler©, Ahrensburg,
 Germany). For each female, we thus calculated an aggregation score, which
 was defined as the total number of pictures in which a female was in the
 chamber next to the group of females for 72 hours. No female was used as
 both a focal and a non-focal individual. We measured the aggregation score
 of 66 females per week for five consecutive weeks (10 died during the
 experiment). Note that the aggregation score (AS) calculated over 24h was
 consistent over the three days of the experiment (Linear models;
 ASday2~ASday1: F<sub>1,310</sub> = 175,52, Adjusted
 R<sup>2</sup> = 0.358, P &lt; 0.001; ASday3~ASday1:
 F<sub>1,310</sub> = 52.85, P &lt; 0.001, Adjusted
 R<sup>2</sup> = 0.132; and ASday3~ASday2:
 F<sub>1,310</sub> = 206.29, P &lt; 0.001, Adjusted
 R<sup>2</sup> = 0. 387) and did not change over the five weeks
 of the experiment (ASday1: F<sub>4,310</sub> = 1.25, P =
 0.292; ASday2: F<sub>4,310</sub> = 0.08, P = 0.990; ASday3:
 F<sub>4,310</sub> = 0.25, P = 0.912; all interactions P
 &gt; 0.141) (Fig S2A). At the end of the 72 hours,
 we isolated each tested female in a Petri dish (5 cm diameter) lined with
 moistened sand for further measurements (see below), while we returned the
 three attractor females to the attraction containers. To avoid
 cross-contamination with biological material between trials
 (<em>e.g.,</em> due to chemical signatures or feces), we
 grounded each plastic arena with a paper sheet that was renewed at each
 trial and cleaned the glass plates with ethanol between each
 trial. For the rest of the experiment, we selected the
 20 females (out of 320) with the lowest aggregation scores (AS from 0 to
 15, n =20; later called low-aggregation females) and the 20 females with
 the highest aggregation scores (AS from 49 to 72, n =20; later called
 high-aggregation females) calculated over 72 hours. These females were
 evenly distributed over the five weeks of the
 experiment, <em>i.e.</em> 4 females per week and category (Fig
 S2B). <em><strong>Behavioral and
 morphological measurements</strong></em> We
 then measured the general activity and boldness of the 40 high- or
 low-aggregation females (Fig S1). At the end of the aggregation test
 (<em>i.e.</em> in the early afternoon), we gently transferred
 each female to a circular arena (diameter 7.8 cm) held between two glass
 sheets and maintained on an infrared light table. We then video-recorded
 females’ activity for 15 minutes in the dark and under infrared light
 (BASLER BCA 1300, Germany; Media Recorder v4.0, Noldus Information
 Systems, Netherlands) – as this species is nocturnal and lucifugous. We
 defined the level of general activity as the total distance walked by a
 female during these 15 minutes [6]. This distance was automatically
 extracted from our videos using the video-tracking software EthoVision XT
 16 (Noldus Information Technology©, Wageningen, Netherlands). We measured
 boldness using the same videos and the same software by extracting the
 time spent by each female in the central area (diameter 5.8 cm) of the
 tested arena. This measurement of boldness is standard in rodents [7, 8]
 and applies to earwigs, as they are highly thigmotactic and typically
 avoid open areas [9]. At the end of the behavioral measurements, we
 analyzed the 72 pictures sampled per female to select the four and four
 females (of the given session/week) with the least and highest aggregation
 scores, respectively. The day after (morning), we
 measured the body condition of the 40 females by first anaesthetizing each
 of them with CO<sub>2</sub> and then measuring their eye
 distance and fresh weight [10]. Using these two measurements, we
 calculated the initial body condition for each female based on the “scaled
 mass index” [11]. In brief, this index standardizes body mass at a fixed
 value of a linear body measurement based on the scaling relationship
 between these measures [12]. Accordingly, this index indicates which mass
 a particular female would have at the average eye distance. We measured
 eye distance to the nearest 0.01 mm using a binocular scope coupled to the
 <em>Leica Application Suite</em> software (Leica Biosystems©,
 Wetzlar, Germany). We weighed each female using an isoCaL Quintix®
 precision balance to the nearest 0.01 mg (Sartorius©, Göttinger, Germany).
 Note that the 8 females from the last week were dissected before weight
 measurement and were thus not included in this measurement.
 <strong>DNA extraction and amplification
 protocol</strong> Immediately after morphological
 measurement (<em>i.e.</em>, less than 24 hours after the last
 aggregation score recording), we extracted the gut of the 40 selected
 females to investigate whether high or low aggregation females had
 different gut microbial communities. The remaining females (and males)
 were used for other experiments not shown in this study. Following the
 protocol detailed in Van Meyel <em>et al.</em> (2021), each
 female was CO<sub>2</sub> anaesthetized and dissected under
 sterile conditions. We immediately flash-freeze each gut with liquid
 nitrogen and stored them individually in Eppendorf tubes at -80°C until
 DNA extraction. Total genomic DNA was extracted using the NucleoMag®
 Tissue extraction kit (Macherey-Nagel™, Düren, Germany) and the V3-V4
 region of the 16S rRNA gene was amplified with the prokaryotic primers
 343F (5′- ACGGRAGGCAGCAG – 3′) and 784R (5′- TACCAGGGTATCTAATC – 3′) [14]
 using the Taq Polymerase Phusion® High-Fidelity PCR Master Mix with GC
 buffer (Qiagen, Hilder, Germany). PCR amplification steps involved an
 initial denaturing step for 30 sec at 98°C, followed by 22 amplification
 cycles (denaturation for 10 sec at 98°C; annealing for 30 sec at 61.5°C;
 extension for 30 sec at 72°C), and ended by a final extension step of 10
 min at 72°C. We checked the amplification success and probe specificity by
 electrophoresis on a 1.5% agarose gel. We amplified all samples in
 duplicates and equally pooled them for a final product of 30 µL. Each step
 involved several blank samples. We verified the inter-sample contamination
 by testing the material used during the dissection and extraction steps
 (<em>i.e.,</em> dissection tools washed in T1 buffer, sample
 tubes empty with T1 buffer). These samples were placed in the same
 extraction plate of the other samples and followed the same PCR protocol.
 Finally, blanks with molecular water were also used during the
 amplification steps. As all blanks were negative, they were not sequenced
 (Figure S3). We sequenced each amplicon by 2x250bp Illumina MiSeq
 technology at the Bio-Environnement platform (University of Perpignan,
 France). <strong>Bioinformatic processing and
 diversity metrics</strong> The obtained sequence
 reads were trimmed and filtered using the quality profiles from the DADA2
 algorithm v1.24.0 [15], cleaned for errors, dereplicated, and inferred
 towards Amplicon Sequence Variants (ASVs) [16]. Chimaeras were removed and
 taxonomy assignment was performed on a count table where forward and
 reverse ASVs were merged and pooled, using the SILVA reference database
 (release 138) [17]. The multiple alignment of the sequences was provided
 with the MAFFT program [18] and we inferred the phylogenetic tree with the
 FastTree 2 tool [19] and Phangorn package v2.8.1 [20]. The table was then
 transformed into a phyloseq object using the phyloseq package v1.40.0 [21]
 of which we removed the organites and eukarya sequences. Finally, we
 normalized all samples at 10 180 sequences (<em>i.e.,</em> the
 minimum sample sum of the dataset) and obtained a final dataset
 constituted by 397 020 sequences (1 593 ASVs) for 20 high- and 19
 low-aggregation females. For each sample, we calculated taxonomic and
 phylogenetic alpha diversity with the Shannon and the Hill1 indices (q =
 1) [22]. We then measured the taxonomic and phylogenetic dissimilarities
 computed on relative abundances of ASVs between pairs of gut microbiota
 using Bray-Curtis and Weighted Unifrac distances, respectively [23]. The
 community homogeneity (<em>i.e.</em> microbial dispersion)
 within low- and high-aggregation status is calculated for each beta
 diversity metric used with the function betadisper from the package
 <em>vegan</em> v2.6.2 [24].
 <strong>Statistical analyses</strong>
 We conducted all statistical analyses using R v4.2.0 [25] loaded
 with <em>car</em> v3.1 [26], <em>lme4</em> v1.1.29
 [27], <em>vegan</em>v2.6.2 [28] and <em>ape</em>
 v5.6.2 [29]. We ran six independent statistical models where the response
 variable was either Shannon (Linear Model, LM), Hill1 (LM), Bray-Curtis
 dispersion (LM), Weighted Unifrac dispersion (LM), Bray-Curtis
 dissimilarity matrix (Permutational Multivariate Analyses of Variances;
 PERMANOVA) or Weighted Unifrac dissimilarity matrix (PERMANOVA). In each
 of these models, we entered the same 5 explanatory variables, namely
 aggregation status (high or low), general activity, boldness, body
 condition (<em>i.e.,</em> scaled body mass index), and
 sampling week. We first included the interaction between sampling week and
 each variable in the different models and then, where appropriate, removed
 it after model simplification by AIC comparison. The PERMANOVA on the
 beta-diversity dissimilarity matrices included 999 permutations and was
 run using the adonis2 function from the package <em>vegan</em>
 v2.6.2 [28]. We performed post-hoc pairwise Adonis tests between sampling
 weeks when it was significant in the PERMANOVA. For each statistical
 model, we checked that all assumptions were fulfilled using the
 <em>DHARMa</em> R package v0.4.6 [24]. There was no
 correlation between the continuous variables entered in each statistical
 model (Pearson tests; P &gt; 0.05 and |R| ≤ 0.28 for all
 pairs). To identify potential microbial discriminants
 of females with high or low-aggregation status, we used a Linear
 discriminant analysis with Effect Size (LEfSe, [30]) based on microbial
 Genera conducted on the Galaxy platform. Finally, to
 test whether high- and low-aggregation females differ in their other
 life-history traits, we performed four t-tests on the general activity,
 boldness, and body condition (<em>i.e.,</em> scaled body mass
 index).
 <strong>References</strong> 
 Wirth T, Le Guellec R, Vancassel M, Veuille M (1998) Molecular
 and reproductive characterization of sibling species in the European
 earwig (<em>Forficula auricularia</em>). Evolution 52:260–265.
 https://doi.org/10.1111/j.1558-5646.1998.tb05160.x González-Miguéns R, Muñoz-Nozal E, Jiménez-Ruiz Y, et al (2020) Speciation patterns in the <em>Forficula auricularia</em> species complex: cryptic and not so cryptic taxa across the western Palaearctic region. Zoological Journal of the Linnean Society 190:788–823. https://doi.org/10.1093/zoolinnean/zlaa070 Meunier J, Wong JWY, Gómez Y, et al (2012) One clutch or two clutches? Fitness correlates of coexisting alternative female life-histories in the European earwig. Evol Ecol 26:669–682. https://doi.org/10.1007/s10682-011-9510-x Van Meyel S, Devers S, Meunier J (2022) Earwig mothers consume the feces of their juveniles during family life. Insect Science 29:595–602. https://doi.org/10.1111/1744-7917.12941 Van Meyel S, Meunier J (2022) Costs and benefits of isolation from siblings during family life in adult earwigs. Animal Behaviour 193:91–99. https://doi.org/10.1016/j.anbehav.2022.09.003 Merleau L-A, Larrigaldie I, Bousquet O, et al (2022) Exposure to pyriproxyfen (juvenile hormone agonist) does not alter maternal care and reproduction in the European earwig. Environ Sci Pollut Res 29:72729–72746. https://doi.org/10.1007/s11356-022-20970-z Archer J (1973) Tests for emotionality in rats and mice: A review. Animal Behaviour 21:205–235. https://doi.org/10.1016/S0003-3472(73)80065-X Stratton JA, Nolte MJ, Payseur BA (2021) Evolution of boldness and exploratory behavior in giant mice from Gough Island. Behav Ecol Sociobiol 75:65. https://doi.org/10.1007/s00265-021-03003-6 Rankin SM, Palmer JO (2009) Dermaptera. In: Encyclopedia of Insects. Elsevier, pp 259–261 Thesing J, Kramer J, Koch LK, Meunier J (2015) Short-term benefits, but transgenerational costs of maternal loss in an insect with facultative maternal care. Proc R Soc B 282:20151617. https://doi.org/10.1098/rspb.2015.1617 Kramer J, Meunier J (2016) Maternal condition determines offspring behavior toward family members in the European earwig. BEHECO 27:494–500. https://doi.org/10.1093/beheco/arv181 Peig J, Green AJ (2010) The paradigm of body condition: a critical reappraisal of current methods based on mass and length: The paradigm of body condition. Functional Ecology 24:1323–1332. https://doi.org/10.1111/j.1365-2435.2010.01751.x Van Meyel S, Devers S, Dupont S, Dedeine F (2021) Alteration of gut microbiota with a broad-spectrum antibiotic does not impair maternal care in the European earwig. Peer Community in Evolutionary Biology Muyzer G, de Waal EC, Uitterlinden AG (1993) Profiling of complex microbial populations by denaturing gradient gel electrophoresis analysis of polymerase chain reaction-amplified genes coding for 16S rRNA. Appl Environ Microbiol 59:695–700. https://doi.org/10.1128/aem.59.3.695-700.1993 Callahan BJ, McMurdie PJ, Rosen MJ, et al (2016) DADA2: High-resolution sample inference from Illumina amplicon data. Nat Methods 13:581–583. https://doi.org/10.1038/nmeth.3869 Glassman SI, Martiny JBH (2018) Broadscale ecological patterns are robust to use of Exact Sequence Variants versus Operational Taxonomic Units. mSphere 3:e00148-18. https://doi.org/10.1128/mSphere.00148-18 Quast C, Pruesse E, Yilmaz P, et al (2012) The SILVA ribosomal RNA gene database project: improved data processing and web-based tools. Nucleic Acids Research 41:D590–D596. https://doi.org/10.1093/nar/gks1219 Katoh K (2002) MAFFT: a novel method for rapid multiple sequence alignment based on fast Fourier transform. Nucleic Acids Research 30:3059–3066. https://doi.org/10.1093/nar/gkf436 Price MN, Dehal PS, Arkin AP (2010) FastTree 2 – Approximately maximum-likelihood trees for large alignments. PLoS ONE 5:e9490. https://doi.org/10.1371/journal.pone.0009490 Schliep KP (2011) phangorn: phylogenetic analysis in R. Bioinformatics 27:592–593. https://doi.org/10.1093/bioinformatics/btq706 McMurdie PJ, Holmes S (2013) phyloseq: an R package for reproducible interactive analysis and graphics of microbiome census data. PLoS ONE 8:e61217. https://doi.org/10.1371/journal.pone.0061217 Chao A, Chiu C-H, Jost L (2010) Phylogenetic diversity measures based on Hill numbers. Phil Trans R Soc B 365:3599–3609. https://doi.org/10.1098/rstb.2010.0272 Lozupone C, Knight R (2005) UniFrac: a new phylogenetic method for comparing microbial communities. Appl Environ Microbiol 71:8228–8235. https://doi.org/10.1128/AEM.71.12.8228-8235.2005 Hartig F (2022) DHARMa: Residual diagnostics for hierarchical (multi-Level / mixed) regression models. R package version 0.4.6 R Core Team (2022) R: A language and environment for statistical computing Fox J, Weisberg S, Fox J (2019) An R companion to applied regression, 3rd ed. Sage, Thousand Oaks, Calif Bates D, Mächler M, Bolker B, Walker S (2015) Fitting linear mixed-effects models using lme4. J Stat Soft 67:. https://doi.org/10.18637/jss.v067.i01 Oksanen J, Simpson G, Blanchet F, et al (2022) vegan: community ecology package Paradis E, Schliep K (2019) ape 5.0: an environment for modern phylogenetics and evolutionary analyses in R. Bioinformatics 35:526–528 Segata N, Izard J, Waldron L, et al (2011) Metagenomic biomarker discovery and explanation. Genome Biol 12:R60. https://doi.org/10.1186/gb-2011-12-6-r60
description: # With or without you: Gut microbiota does not predict aggregation
 behavior in European earwig females
 [https://doi.org/10.5061/dryad.z8w9ghxmm](https://doi.org/10.5061/dryad.z8w9ghxmm) This repository provides all data necessary to reproduce the study of Cheutin MC, Lelclerc B, Meunier J. 2024 from *With or without you: Gut microbiota does not predict aggregation behavior in European earwig females*. Behavioral Ecology. Briefly, it contains the following data **(1)** a zip folder ("**Aggregation_reads.gz.zip"**) that contains the original 16s rRNA libraries (R1R2) not processed .fastq.gz format; **(2)** A sample table containing all information related to the female earwig hosts (behavioral and biological traits) (**"env_data.txt"**), **(3)** A metadata table containing each aggregation score per hour per individual (**"meta_social.txt"**), **(4)** an R phyloseq object (**ps_rff.RData**) used for the analyses. We also provide two R scripts: **(1)** for *sequence processing* (called **"Setup.R"**) to proceed with the reads from the **Aggregation_reads.gz.zip** folder and create the *rarefied phyloseq* file (**ps_rff.RData**) used for further analyses; **(2)** A R script for further *data analyses* (**"Analyses.R"**) that is divided in two: **(2A)** *Aggregation analyses* that allows to calculate the aggregation scores (AS) for the 72h, test its repeatability and its distribution through the earwig female population. This part calls for the table (**"meta_social.txt"**). **(2B)** *Microbial analyses* that allow the analysis of the microbial composition, measure the alpha and beta diversities, test their determinants, and analyse the bacterial biomarkers associated with aggregation profiles (i.e., bacteria that are indicators of low- versus high-aggregative earwig hosts). This part calls for the table (**"env_data.txt"**). ## Description of the data and file structure ### 1. Aggregation\_reads.gz.zip The folder contains the original 16s rRNA libraries (R1R2) not processed in .fastq.gz format. You can unzip these files using the Unarchiver tool. ### 2. env\_data.txt * *ID* : Sample identification of the 39 earwig hosts (categorical) * *sampling_date* : date of the experimental sampling (categorical ; Run 1, Run 2 , Run 3, Run 4 and Run 5) * *D1* : aggregation score at Day 1 (discrete ; counts [0:24]) * *D2* : aggregation score at Day 2 (discrete ; counts [0:24]) * *D3* : aggregation score at Day 3 (discrete ; counts [0:24]) * *distance_moved* : distance moved (cm) by the sample for 15 minutes (continuous) * *cum_dur_center* : cumulative duration (sec) in the center of the activity arena (continuous) * *status* : aggregation status (binary ; Low- or High-aggregative). * total_AS : the aggregation score at the end of the 72H (i.e., sum of D1, D2 and D3) of the experimental measurement (discrete ; counts [0:72]). * fresh_mass : measured fresh mass (mg) (continuous) * eyes_distance : distance between the eyes of the samples (mm) (continuous). * *SMI* : scaled mass index, a proxy for the body mass index. This index standardizes body mass at a fixed value of a linear body measurement (eyes distance) based on the scaling relationship between these measures (continuous). * *R_length* : length (mm) of the right forceps (continuous) * *L_length* : length (mm) of the left forceps (continuous) * *Average_length* : Average length (mm) of the right and left forceps (continuous) ### 3. meta\_social.txt * *ID* : Sample identification of the 320 earwig samples measured for aggregation (categorical) * *sampling_date* : date of the experimental sampling (categorical ; Run 1, Run 2 , Run 3, Run 4 and Run 5) * *attraction* : the value of the attraction related to the position of the earwig observed (binary ; 0: out of the aggregative chamber ; 1 : in the aggregative chamber). * *time_exp* : time of the experimental sampling(discrete ; from 1 to 72) * *day* : day of the experimental sampling (categorical ; D1, D2 or D3) ### 4. ps\_rff.RData This R object is the phyloseq object that contains the following structure : * *otu_table()* : 1593 taxa (ASVs) and 39 samples (count table) * *sample_data()* : 39 samples by 35 sample variables (env_data) * *tax_table()* : 1593 taxa (ASVs) by 8 taxonomic ranks (descriptive taxonomical table) * *phy_tree()* : 1593 tips and 1592 internal nodes (phylogenetic tree of the ASVs) *refseq()* : 1593 reference sequences This file can be loaded with the following command in R v4.2.2 *&gt;load("ps_rff.RData")* but see the below scripts. ## Code/Software R is required to run the following script. Note that the package *phyloseq* might need a recent version of R (here v4.2.2). #### 1. *Setup.R* script The script contains all information regarding the bioinformatic processing in R from the preparation of the working space and the library packaging to the sequence processing in order to create the rarefied phyloseq object ps_rff.RData. This part calls for the R1R2 libraries contained in the folder **Aggregation_reads.gz.zip** and the **env_data.txt** file to create the phyloseq object. #### 2. *Analyses.R* script All analyses are made in R with the **"meta_social.txt"** for the *Aggregation analyses* and the **"ps_rff.RData"** file for the *Microbial analyses*. * I - Aggregation analyses: The code allows the description of the frequency of the aggregation score observed during 72h in the entire population and allows testing its repeatability during the 3 days. * II - Microbial analyses: The code allows to describe graphically the composition of the microbiota between samples at different phylogenetic levels (i.e., microbial Phylum, Family, and Genus). A second part is dedicated to the measure of the alpha and beta diversity. A third part allows the testing of the relation between the microbial diversity of female earwigs and their aggregation status (Low- or High-aggregative hosts) and between microbiota and morphological measurements.