Gaps in our understanding of bacterial species distributions, functions and communities are of particular interest within wildlife ecology because pathogens (bacterial or otherwise) have the potential to regulate host populations [Reference Anderson and May1], which can in turn drive ecosystem dynamics. In addition, zoonoses pose a significant human health risk in much of the world, accounting for 60% of emerging infectious diseases documented since 1940, 70% of which originated in wildlife [Reference Jones2].

Both of these concerns intersect in rodents. Several species of vole and lemming in Europe and North America, including the fossorial water vole Arvicola terrestris in France, exhibit multi-annual population cycles characterized by explosive population growth followed by rapid catastrophic population declines [Reference Saucy3]. These cycles are not completely understood, and pathogens are one potential driver of these cycles that has yet to be fully explored. Rodents can also be reservoir species for zoonotic bacteria around the world [Reference Meerburg, Singleton and Kijlstra4].

The development of high-throughput sequencing and 16s rRNA amplicon sequencing approaches, which allows for the simultaneous identification of multiple bacterial species in large sample sizes, now allows for the study of bacterial communities within hosts and host populations [Reference Razzauti5, Reference Galan6], but methodological questions still remain. For example, tissue selection within hosts is one area that has received little attention. Comparing multiple organ tissues within the same rodent has been done before for single bacterial species [Reference Gundi7, Reference Rossow8], but to our knowledge, the present study is the first to survey multiple organs for their bacterial assemblages.

Our objectives were the following: first, describe the bacterial pathogens present in 13 A. terrestris heart, liver, lung, kidney and spleen tissues. Currently little is known about the bacteria present in this species in the Franche-Comté region of eastern France. Our second objective was to compare bacterial assemblages among organs. We hypothesized that bacterial assemblages would differ in richness and membership among tissue types because the five organs we have targeted differ in both function and potential exposure routes. We were also interested in determining how many organ types are required to obtain complete coverage of the total bacterial assemblage of the five organs. We hypothesized that there would not be complete species turnover between organs, and that a subset of organs would contain the majority of bacterial pathogens found in the pooled assemblage.

We conducted sampling in the communes of Arc-sous-Montenot (46·925° N, 6·005° E) and Censeau (46·802° N, 6·065° E) in the Doubs and Jura departments of eastern France in November of 2014. The two sampling locations are approximately 14 km apart, and separated by a mix of pasture, cultivated agricultural land and forest. Population abundance at the two sites was assessed by at least 1·5 km of 5 m-wide transects [Reference Giraudoux9]; the proportion of 10 m intervals within the transects positive for recent vole activity (tumuli) were 91% and 7% in Arc-sous-Montenot and Censeau, respectively. We collected nine adult males from Arc-sous-Montenot and four adult males from Censeau as part of a larger sampling program (metadata available through the dat@osu platform [10]). We trapped the animals using unbaited Sherman live traps spaced at least 10 m apart, and euthanized them by cervical dislocation. The 13 animals’ weight ranged from 75 g and 97 g (mean 84·3 g), and none of the animals were reproductively active at time of collection (testes were retracted into the abdomen). Portions of the heart, lungs, liver, spleen and kidneys were stored individually in 2 ml tubes of RNAlater (Thermo Fisher Scientific, Illkirch-Graffenstaden France) and subsequently frozen at −20 °C. To prevent cross-contamination between animals during dissection, we systematically alternated the use of two sets of dissecting instruments, which were decontaminated using bleach, water and 95% ethanol between animals [Reference Auffray11]. We did not decontaminate dissecting instruments between organs within animals, and therefore did not eliminate the possibility of cross-contamination between organs within animals. Animals were treated in accordance with European Union guidelines and legislation (Directive 86/609/EEC). The Chrono-environment laboratory received approval from the Comité d'Ethique Bisontin en Expérimentation Animale (CEBEA N°58) for the sampling of rodents and the storage and use of their tissues. The rodent species investigated in this study does not have protected status (see IUCN and CITES lists), and is listed as a pest, subject to control, under Article L201-1 of the Code Rural et de la Pêche Maritime of French law.

We extracted DNA from 65 tissue samples (13 individuals × 1 tissue sample for each of the five organs collected per animal) using the DNeasy® 96 tissue kit (Qiagen, Hilden, Germany) following the manufacturer instructions. We measured DNA concentration for each sample using a NanoDrop 8000 UV-Cis spectrophotometer (Thermo Fisher Scientific, France); all samples had a DNA concentration higher than 10 ng/μl.

We sequenced a 251 bp fragment of the V4 region of the 16S rRNA gene following Kozich et al.’s procedure [Reference Kozich12] with modifications detailed by Galan et al. [Reference Galan6] for PCR amplification, indexing, multiplexing and demultiplexing and taxonomic identification using the SILVA SSU database as a reference (see Supplementary Material 1, available on the Cambridge Core website, for details). Sequences were grouped into Operational Taxonomic Units (OTUs) using a 97% similarity threshold, and only those OTUs identified as potentially pathogenic were included in the subsequent analysis (see Supplementary Material 1 for details). In the rest of this article, the terms ‘OTU’ and ‘bacteria’ are used to describe potentially pathogenic agents only.

All statistics were conducted using R (version 3.1·2, R Foundation for Statistical Computing, Vienna, Austria). We calculated the average number of OTUs detected in each tissue sample (OTU richness) as well as the average OTU richness in each animal when tissues were pooled together. We used the glmer function in the lme4 package [Reference Bates13] to fit a generalized linear mixed model to determine if OTU richness differs between organs. We used a Poisson distribution with log link function, animals nested in locations were included as random factors, and the model was fit using maximum likelihood (Laplace approximation). To determine if pooling organ data significantly increased OTU richness, we used bootstrapping (999 iterations) to estimate the average difference (with 95% confidence intervals) between the heart and pooled OTU richness within animals.

We calculated Jaccard's presence/absence dissimilarity indices between each organ within animals to quantify bacterial assemblage differences, and we used the adonis function in the vegan package [Reference Oksanen14] to perform a PERMANOVA (permutational multivariate analysis of variance) on the Jaccard dissimilarities to determine if bacterial assemblages within organs were non-random. We stratified organs by animal to account for any between-animal variation, and included a location term and interaction term between location and organ as well. We used the betadisper function in the vegan package to check for homogeneity of dispersions of the groups, and set the number of permutations to 999. We also performed a PERMANOVA analysis to determine if bacterial assemblages differ between animals (see Supplementary Material 4 for details).

We performed an OTU accumulation analysis to assess the effect of including more than one organ in calculating animal OTU richness. Organs were included in the OTU richness calculation for each animal in order of decreasing average OTU richness (first the heart, then the kidneys, then the spleen, then the liver, then the lungs). Means and 95% confidence intervals were estimated via bootstrapping (999 iterations).

A total of 24 OTUs representing at least 17 genera were included in the final analysis (Table 1). One OTU could not be assigned to a genus, its family classification (Pasteurellaceae) is listed instead.

Table 1. OTU genera included in the analysis, and their prevalences (expressed as a percentage in all samples and each organ)

One OTU could not be classified to genus level, and is instead designated by family.

^a Two different OTUs.

b Two different OTUs.

Prevalence of individual OTUs ranged from 2% to 37% of all 65 samples, with prevalence within an organ ranging from 0% to 54% (Table 1). Average OTU richness within each organ ranged from 1·5 ± 0·4 (mean ± standard error) to 2·5 ± 0·4 OTUs/organ (Table 2) and in general did not differ significantly between organs with the exception of the average lung OTU richness, which is marginally significantly lower than the average heart OTU richness (Table 3). Location did not account for any variance in the model, and was discarded. When pooling organs together within animals, average OTU richness was 4·7 ± 0·6 OTUs, which is 2·15 ± 0·84 OTUs (bootstrapped mean and 95% confidence interval) greater than the average OTU richness of the heart.

Table 2. Average OTU richness and standard error for each organ, and for organs pooled within animals

Table 3. Parameter estimates (log scale) and standard errors for fixed effects and variance and standard deviation for random effects of a generalized linear mixed effects model for average OTU richness within organs

The model was fit using maximum likelihood, a Poisson distribution and log link function.

We determined that bacterial assemblages differ between organs (PERMANOVA, pseudo-F _4,51 = 1·37, P = 0·001) and between locations (pseudo-F _1,51 = 11·77, P = 0·003). The organ–location interaction term was also significant (pseudo-F _4,51 = 1·29, P = 0·027), suggesting that location effects are not the same across organs (see Supplementary Material 3 for additional analysis).

Adding organs increased the average OTU richness of each animal (Fig. 1); adding the kidney increased the average animal-level OTU richness by one OTU, and the addition of the spleen increased the average OTU richness by 0·8 OTU. Including the least-rich organs last, the liver and lung, still increased average OTU richness by 0·2 OTUs each.

Fig. 1. Bootstrapped operational taxonomic unit (OTU) accumulation curve of the number of potentially pathogenic bacterial OTUs detected within 13 male Arvicola terrestris with increasing number of organs screened. Organs are ordered by average OTU richness, with the heart having the highest average OTU richness, and the lung having the lowest. The grey band is bootstrapped 95% confidence intervals (999 iterations). Animals were collected in Franche-Comté, France, in November 2014.

We failed to detect differences in average OTU richness between organs, but we determined that bacterial assemblages vary between organs within individuals; if all of the organs within an animal had the same bacterial assemblages, we would expect the OTU accumulation analysis to show no increase in average OTU richness with increasing number of organs. This result is further supported by our PERMANOVA results, as differences in assemblages do appear to be mediated by organ type. This suggests that a holistic approach is better suited to surveying bacteria in this species, as bacterial assemblages within different organs within the same animal vary in both richness and membership. These results also suggest that cross-contamination between organs during the dissection and collection process was not significant; if significant cross-contamination had occurred, we would expect bacterial assemblages obtained from different organs within an animal to be similar.

The differences we detected in bacterial assemblages within organs may be due in part by spatial and temporal dynamics of bacterial infections. Deng et al. [Reference Deng15] report that while Bartonella infections result in high Bartonella densities in the spleen, Bartonella in the liver is detectable only early in the infection and for a brief period of time, and similar infection dynamics have been observed for Borrelia burgdorferi in mammals [Reference Barthold16]. In addition, infectious bacteria may target organ tissues but the spatial and numerical distribution of bacteria within host cells can vary with host. For example, Salmonella infections in mice can be characterized by either the widespread but low-abundant presence of bacteria within host cells [Reference Sheppard17], or a more localized infection with few host cells harboring many bacteria [Reference Salcedo18]. This potential for different infection patterns at the cellular level may be of particular importance here as we did not homogenize our organs before taking tissue samples; the apparent increase in OTU richness with increasing number of organs may in fact be a consequence of resampling the same animal five times, rather than differences in assemblages within each tissue type.

Host immuno-competence, body condition and genetics can also contribute to infection dynamics [Reference Beldomenico and Begon19, Reference Turner20] and may drive some of the variation we observed between animals from the same population.

The conclusions presented here are of course limited by the small sample size; OTU accumulation analyses for individual organs (Supplementary Material 2) strongly suggest that if we were to include additional animals in our analysis, we would detect more OTUs. Additional tissues like lymph nodes, the brain, the digestive tract and skin lesions can harbor bacteria that may also be of interest, and experimenting with tissue homogenization and sub-sampling could provide insight into the results we obtained here. Determining if blood would yield results similar to the heart or other organs would also be of interest, given that all of the tissues we considered here are permeated by blood.

Our results demonstrate that the number of organs sampled influences the power to detect trends in bacterial pathogen assemblage composition. These results can inform sampling decisions in public health and wildlife ecology.