Whether, when, and how genetic diversity buffers individuals and populations against infectious disease risk is a critical and open question for understanding wildlife disease and zoonotic disease risk. Several, but not all, studies have found negative relationships between infection and heterozygosity in wildlife. Since they can host multiple zoonotic infections, we sampled a population of wild deer mice (Peromyscus maniculatus), sequenced their genomes, and examined their fecal samples for coccidia and nematode eggs. We analyzed coccidia infection status, abundance, and coinfection status in relation to per-locus and per-individual measures of heterozygosity, as well as identified SNPs associated with infection status. Since heterozygosity might affect host condition, and condition is known to affect immunity, it was included as a co-variate in the per-individual analyses and as response variable in relation to heterozygosity. Not only did coccidia-infected individuals have lower levels of genome-wide per-locus diversity across all metrics, but we found an inverse relationship between genomic diversity and severity of coccidia infection. We also found weaker evidence that coinfected individuals had lower levels of private allelic variation than all other groups. In the per-individual analyses, relationships between heterozygosity and infection were marginal but followed the same negative trends. Condition was negatively correlated with infection, but was not associated with heterozygosity, suggesting that effects of heterozygosity on infection were not mediated by host condition in this system. Association tests identified multiple loci involved in the inflammatory response, with a particular role for NF-κB signaling, supporting previous work on the genetic basis of coccidia resistance. Taken together, we find that increased genome-wide neutral diversity, the presence of specific genetic variants, and improved condition positively impact infection status. Our results underscore the importance of considering host genomic variation as a buffer against infection, especially in systems that can harbor zoonotic diseases.

Field sampling methods

Deer mice (Peromyscus maniculatus) were live-trapped at the Bernard Field Station in Claremont, CA from June 27 to August 9, 2019. Eight 40 × 40 m grids with traps spaced 10 m apart (25 traps/grid) were established in paired locations throughout the sage-scrub habitat. At each grid, Sherman traps were baited with seeds for three nights in a row. Three-night trapping sessions were repeated at each grid every two weeks. Upon their first capture each session, mice were weighed, measured, and a retro-orbital blood sample was collected. Each mouse was also given an individually-numbered ear tag upon its first capture. Females were considered reproductive if they were pregnant or lactating, and males were classified as reproductive based on testes size. Fecal samples from each individual were collected from the traps, weighed, and stored in neutral-buffered formalin. Blood samples were kept on ice, then spun at 2000 g for 10 min to separate plasma. Plasma and pellet fractions were frozen separately at -80°C. Host trait and parasite infection data were used from each individual’s initial capture, while the blood samples sequenced sometimes came from subsequent recaptures. 

Parasite and condition measurements

To determine body condition, we regressed weight with body length then used the residuals as a measure of condition. The number of coccidia oocysts and helminth eggs per g of feces were counted using a modified McMaster egg fecal counting protocol (MAFF 1980). There were 9 uninfected and 22 coccidia-infected individuals. Following the methods of Ezenwa et al. (2021), we used a model-based approach to identify statistically significant clusters in the oocyst intensity data using the mclust package (Fraley et al. 2021). We then classified infected individuals into two groups based on these clusters [Low (<100 oocysts per gram feces (opg)): n = 11; High (>100 opg): n = 10]. For two individuals, feces weight data were not collected, so only infection status, not abundance, could be determined. One of these two individuals only had a single coccidia oocyst detected, so it was included in the Low infection severity category. Lastly, we classified individuals as whether or not they were coinfected with coccidia (Eimeria sp.) and gastrointestinal nematodes (Aspicularis americana (a pinworm) and 3 stryongylid (hookworm) species) (Uninfected: n = 9; Low coccidia (<100 opg): n = 9; High coccidia (>100 opg): n = 7; Coinfected: n = 5). 

Sequencing and bioinformatics workflow

DNA was extracted from heparin-treated blood pellet samples with the Qiagen DNEasy Blood and Tissue Extraction Kit. DNA was quantified on a Qubit™ 2.0 fluorometer using a dsDNA High Sensitivity Assay Kit (Invitrogen). Sequencing libraries were prepared by Admera Health with a KAPA Hyper Prep Kit. Individually barcoded libraries were sequenced on an Illumina NovaSeq S4 (paired end, 2 × 150 bp), at a target depth of ~200 million PE reads per sample.

We conducted an initial assessment of read quality using FASTQC v0.11.8 with default values (Andrews 2010). TRIMMOMATIC v0.39 was then used to trim adapter sequences, the leading and trailing three bases, and low quality sequences ( leading:3 trailing:3 slidingwindow:4:15 minlen:36; Bolger et al. 2014). Filtered sequences were mapped to the P. maniculatus genome (GCF_000500345.1_Pman_1.0_genomic.fna) using BWA mem v0.7.17 with default values (Li). After sorting (“sort”) and indexing (“index”) BAM files with SAMTOOLS v1.14 (Li et al. 2009), PICARD v 2.26.1 was used to remove duplicates (“MarkDuplicates”) and add read groups ("AddOrReplaceReadGroups"; validation_stringency = lenient; Broad Institute 2019).

We called variants in the parallel mode of FREEBAYES v1.3.5 with 16 threads, resulting in an initial total of 153,221,207 variant sites (Garrison and Marth 2012). To retain only high quality SNPs, we first used vcftools 0.1.17 to filter variant sites with > 5% missing data (–max-missing 0.95), quality scores < 40 (–minQ 40), minimum depth <8 (–minDP 8), and a minor allele frequency of at least 1.5% (–maf 0.015) (Danecek et al. 2011). We next filtered for allele balance, or the proportion of reads supporting a variant, to reduce paralogs and errors using the vcffilter tool in the vcflib package v1.0.1 (Garrison et al. 2021). We removed sites with an allele balance <0.25 or >0.75, with the exception that we retained allele balances close to zero to keep variants where all individuals were homozygous for one or the other allele. Vcftools v0.1.17 was then used to retain only biallelic SNPs with a genotype quality >20 (–min-alleles 2, –max-alleles 2, –minGQ 20, –remove-indels) and remove any loci that were out of Hardy–Weinberg Equilibrium (P<0.01; –hwe 0.01) in Uninfected individuals (Danecek et al. 2011;). Linkage disequilibrium filtering was performed by pruning out sites with r2>0.6 in 1000 site windows with bcftools v1.11 (+ prune -m 0.6 -w 1000; Li 2011; Danecek et al. 2021). After filtering, 573,838 high quality SNPs were retained for downstream analyses.

Genetic diversity analyses

Per-individual genetic diversity metrics were calculated in vcftools v0.1.17 and included HO, HE, and FIS. Per-individual measures of diversity average across all loci in the genome for each individual (n = 31 individuals).

BFS 2019 - rodent infection heterozygosity data.csv contains capture data, infection data and heterozygosity data.

Capture data: Trapping grid (GRID), grid site (Grid2), trap number (Trap), Species (PM = Peromyscus maniculatus), Sex (F = female, M = male),  Age (A = adult, SA = subadult), reproductive characteristics (VAG = perforate vagina, NIP = enlarged nipples indicative of lactation, PREG = pregnant, TES = enlarged testes, Repro = reproductive status, total length in mm, body length in mm, and weight in g. 

Infection data: Coccidia oocysts per gram of feces at initial capture (Init.Coccidia.EPG, and log transformed in logCoccidiaEPG), nematode eggs per gram of feces at initial capture (Init.Nem.EPG and log transformed in logNemEPG), binary infected with coccidia or not (Infection.Status), coinfection status (Coinf.status), Coccidia severity (Coccidia.group), and length-weight residual condition (condition)

Individual-level heterozygosity data: Number of sites analyzed (N_SITES), expected heterozygosity (He), standardized observed heterozygosity (Ho_standard), and standardized inbreeding coefficient (Fis_standard)

Pman_13Jul21.vcf.gz

This file contains SNP data.