Investigation of Genetic Susceptibility to Blastomycosis Reveals Interleukin-6 as a Potential Susceptibility Locus

Blastomycosis is a potentially life-threatening infection caused by the fungus Blastomyces dermatitidis. As with related fungal diseases, blastomycosis is noted to affect some populations more than others. These patterns of illness are often not related to predisposing conditions or exposure risks; thus, genetic differences are thought to underlie these health disparities. People of Hmong ancestry in Wisconsin are at elevated risk of blastomycosis compared to the general population. We studied the genetic codes of Hmong blastomycosis patients and identified candidate sites in their genomes that may explain their susceptibility to this infection. We further studied one particular region of the genome that is involved with the immune processes that fight B. dermatitidis. Our work revealed population differences in the response to fungi. A better understanding of the genetic underpinnings of susceptibility to infectious diseases has broader implications for community health, especially in the paradigm of personalized medicine.

genetic polymorphisms associated with susceptibility to dimorphic fungal infection have been reported, and variants that confer risk of infection are understood to be major mutations of Mendelian loci. Reports of such variants account for only a few individual patient's cases (3,4).
The few genetic studies of Hmong populations have revealed that they are highly distinct from those of other East Asian populations. Overall homozygosity is high in the Hmong (5,6), consistent with historically high rates of within-group marriage and recent population bottlenecks due to repeated migrations, geographic dispersal, and population losses from military conflict (7). This history is expected to lead to increased levels of autozygosity, that is, homozygosity in which two alleles are identical by descent. We hypothesize that homozygous regulatory variant(s) common to the WI Hmong affect their immune responses to fungi, leading to an elevated incidence of blastomycosis.
In this report, we used a bioinformatic strategy to identify polymorphisms that may account for the increased risk of blastomycosis among the WI Hmong and functionally validated our findings with in vitro studies of human cells and a murine model of blastomycosis. Briefly, we used a disease mapping strategy designed to exploit the high levels of autozygosity in the study population, as has been used previously to identify homozygous alleles associated with various diseases (8). We then leveraged human genetic resources characterizing variants across the genome to prioritize variants that may mediate susceptibility to infection. We uncovered a set of 25 polymorphisms near the gene encoding interleukin-6 (IL-6), IL6. We focused on these variants in part due to the role of IL-6 in nurturing IL-17 responses, which are protective in murine models of blastomycosis (9). We report significant differences in cytokine production between populations, with Hmong donors producing less IL-6 than do their European counterparts and also less IL-17 by memory CD4 ϩ T cells in response to fungi. Finally, we formally demonstrate the importance of IL-6 in the development of adaptive immunity and resistance in a murine model of pulmonary infection with Blastomyces dermatitidis.

RESULTS AND DISCUSSION
Analysis of Hmong blastomycosis patient genomes reveals 113 candidate susceptibility variants. To identify variants that may affect susceptibility to infection, we performed whole-genome sequencing (WGS) on 11 Hmong volunteers. The mean coverage depths across the genome for each of the 11 samples ranged from 11 to 32ϫ, with the mean being 22ϫ, and 9,825,483 variant sites were identified against the reference genome (Fig. 1A). One donor, the relative of a blastomycosis patient, was flagged from further analysis because they did not have a confirmed history of blastomycosis. Of the 10 remaining donors with a history of confirmed blastomycosis, one additional donor was flagged because their risk for blastomycosis may have been modified by comorbidity (medical history of organ transplantation). Thus, WGS data from nine Hmong donors with a history of treatment for blastomycosis were used in the subsequent disease mapping strategy.
Operating under the hypothesis that susceptibility loci would be enriched in regions of autozygosity, we mapped runs of homozygosity (ROH). Consistent with other worldwide populations (10), the majority of ROH identified in Hmong cases are small (Ͻ49,277.2 bp using 3 clusters; see Fig. S1 in the supplemental material). We identified 27,514 consensus ROH (defined as regions where at least seven of the nine donors have overlapping ROH segments), which comprised 638,091,042 bp of the autosomes (Fig. 1A). The median consensus ROH was 8,773 bp in length, with the largest identified being 11,100,866 bp. ROH segments within 500 bp of another were merged prior to variant annotation and filtering (1,176 segments). We found 2,010,254 variants within these regions, of which 402,206 were fixed among all nine Hmong donors (Fig. 1A).
These fixed single nucleotide polymorphisms (SNPs) were filtered based on frequencies of the alleles in European populations, as Hmong-Americans exhibit higher incidence rates of blastomycosis than do European-Americans (1). We imposed a threshold of 20%, resulting in 6,037 SNPs (Fig. 1A). We used three criteria to prioritize variants with a potential influence on fungal susceptibility ( Fig. 1A; see also Materials and Methods) and identified 63 consensus ROH that harbor candidate susceptibility variants (Fig. 1B). Of 113 candidate susceptibility variants, 59 are annotated as regulatory by the Ensembl Variant Effector Predictor, 56 are within 200 kb of the transcription start site of an antifungal immunity gene, and 38 are demonstrated expression quantitative trait loci (eQTL) in the GTEx project ( Fig. 1A to C and Table S1).
IL6 and AS-IL6, a long noncoding RNA, represent a candidate susceptibility locus. We manually inspected candidate variants and performed functional validation for one of the most promising loci. The gene encoding the cytokine IL-6 (IL6) is   Fig. 1B and 2A). With 25 SNPs within 200 kb of the transcription start site, IL6 is the gene associated with the most candidate susceptibility variants in our analysis (Table S1). The consensus ROH overlapping the IL6 gene harbors 24 candidate susceptibility variants ( Fig. 2A); the vast majority of consensus ROH do not harbor candidate variants, and no other consensus ROH harbors more than 4 candidate variants (median, 1 variant; Fig. 1D). Of note, 22 of the 25 SNPs identified near IL6 are associated with expression of a long noncoding RNA (lncRNA) that overlaps IL6 on the antisense strand (AS-IL6; AC073072.5 in the GTEx project). This likely reflects linkage disequilibrium between these SNPs (Fig. 2B). Nevertheless, AS-IL6 is the gene associated with the most candidate susceptibility variants identified by our two unbiased criteria (Table S1). In addition to these features at the IL6 locus, the encoded cytokine participates in the development of IL-17 responses, which are essential for host resistance in murine models of blastomycosis (9), lending biological plausibility to a role for the IL6 locus in patients.
To further investigate the frequency of candidate SNPs at the IL6 locus, we sequenced a subset of SNPs from a second cohort of Hmong donors (n ϭ 10) that did not report a history of blastomycosis. Though these SNPs are not fixed in the WI Hmong at large, they are at a higher frequency in the Hmong than in the East Asian population as a whole (allele frequency in Hmong [AF Hmong ], 0.975; allele frequency in European ancestry [AF EAS ], 0.783; Fig. 2C). This haplotype is rare in European populations (AF EUR ϭ 0.048; Fig. 2C).
AS-IL6 is reported to play an important role in facilitating IL-6 production via H3K27ac in the IL6 promoter (11). Intriguingly, one candidate SNP (rs1800796) is transcribed in AS-IL6, and an additional two candidate SNPs (rs1524107 and rs2066992) are upstream of AS-IL6 in predicted transcription factor binding sites (SNP FuncPred, https://snpinfo.niehs.nih.gov/snpinfo/snpfunc.html). All three of these SNPs are reported eQTLs for AS-IL6 but not for IL6 in the GTEx project; a high degree of variability in IL6 expression in postmortem tissues may interfere with the detection of eQTLs for this gene. Based on these data, we hypothesize that candidate susceptibility variants we have identified regulate IL6 expression, possibly indirectly by impacting the regulatory lncRNA AS-IL6. In support of this hypothesis, we found a strong correlation between the AS-IL6 expression and IL6 when we reanalyzed published transcriptomic data from monocytes stimulated with lipopolysaccharide (LPS) (R 2 ϭ 0.91, P ϭ 0.002) (12), which may reflect our hypothesis or a shared promoter (Fig. 2D).
Functional validation of ethnic difference in cytokine production. While IL-17 responses have been implicated in murine models of B. dermatitidis infection, they have not been studied in humans with blastomycosis. Recognizing a unique opportunity to study the role of IL-6 and related IL-17 responses in human populations with known susceptibility, we quantified cytokine responses by Hmong and European donors' cells. We were not able to reach Hmong blastomycosis patients for re-collection of primary cells; to circumvent this limitation, we created immortalized B-lymphoblastoid cells (B-LCLs) from Hmong and European subjects. B-LCLs from Hmong donors produced significantly less IL-6 than did B-LCLs from European donors in both the unstimulated (P ϭ 0.002) and stimulated (P Ͻ 0.001) states (Fig. 3A). We also quantified IL-8 production as a control to see if there is a generalized difference in cytokine production; in contrast to IL-6, Hmong subject cells produced significantly more IL-8 upon stimulation than did European subject cells (Fig. 3A). Differences in IL-6 production remained significant (P Ͻ 0.05) regardless of whether the two donors (one Hmong and one European) that are heterozygotic at SNPs near IL6 were included or excluded from analysis. The difference in IL-8 production upon stimulation was similar between the two donor groups if the two heterozygotic donors are excluded (P ϭ 0.052). B-LCLs failed to produce cytokine when stimulated with fungal antigens (data not shown).  Genetic Susceptibility to Blastomycosis ® Our findings of an ethnic difference in IL-6 production, when correlated with genotype at the IL6 locus, fit with past reports of human genetic differences in IL-6 production. Indeed, one candidate SNP (rs1800796) has been repeatedly associated with serum levels of IL-6 (13-16). These studies involve donors from different ethnicities and report contradicting relationships between genotype and IL-6 levels, perhaps indicating epistatic effects of these SNPs.
We also sought to measure antifungal immune responses in our donors downstream of IL-6. We surveyed memory CD4 ϩ T cell responses to Candida albicans because it is likely that most healthy adults have been exposed to this fungus, and we hypothesized that deficits in IL-6 responses to fungi would affect this response. Hmong donor cells produced significantly less IL-17 than did European donor cells in response to C. albicans (P ϭ 0.033; Fig. 3B); Hmong donors also produced less IL-17 in the unstimulated state (P ϭ 0.012). Hmong donor cells likewise had lower expression of RAR-related orphan receptor gamma t (ROR␥t), the hallmark transcription factor of IL-17-producing T cells, in response to C. albicans than did European donor cells (Fig. 3C  and D). The donor groups did not differ in IL-17 levels or ROR␥t expression in response to either direct T cell receptor (TCR) ligation (␣-CD3, ␣-CD28, and phytohemagglutinin [PHA]) or the vaccine antigen tetanus toxoid (Fig. 3B to D). These data suggest that impaired IL-6 production in the Hmong subjects may affect IL-17-producing T cells in response to fungi more generally. We did not detect differences in IFN-␥ responses to heat-killed Candida or TCR ligation but did find that cells from Hmong donors produced less gamma interferon (IFN-␥) than did European donors in response to tetanus toxoid (Fig. 3B). Our genomic analysis of Hmong blastomycosis patients leads us to hypothesize that these IL-17 differences by memory CD4 ϩ T cells in response to fungal antigen are downstream of IL-6. To investigate whether ethnic differences in innate (myeloid) IL-6 or adaptive (IL-6 dependent) IL-17-producing T cell responses would be more likely to account for patient susceptibility to blastomycosis, we studied the loss of IL-6 signaling in murine models of innate and adaptive immunity to B. dermatitidis. IL-6 is required for the development of protective T cells and neutrophil responses in a murine model of blastomycosis. IL-6 promotes the development of T helper 17 (Th17) responses (17), and murine models of antifungal immunity have demonstrated that Th17 cells foster acquired resistance to experimental pulmonary blastomycosis (9). We investigated the contribution(s) of IL-6 to innate and acquired resistance in the murine model of blastomycosis.
We found that IL-6 had a requisite role in the development of adaptive immunity and resistance to B. dermatitidis. In a vaccine model, IL6 Ϫ/Ϫ mice failed to acquire resistance to lethal experimental challenge, as evidenced by significantly higher fungal burden than in wild-type (WT) mice (Fig. 4A). IL6 Ϫ/Ϫ mice recruited significantly fewer total CD4 ϩ T cells and IL-17-producing CD4 ϩ T cells to the lung upon fungal challenge than did vaccinated WT mice (Fig. 4B). IFN-␥ responses were similar between the groups (Fig. S2A). Splenocytes from the IL6 Ϫ/Ϫ mice likewise revealed a trend toward less IL-17 than those from WT mice upon ex vivo stimulation with B. dermatitidis antigen (Fig. 4B). Upon excluding one outlier in the WT vaccinated group (Grubbs' test, P Ͻ 0.05), this difference in IL-17 production between groups reached significance (P ϭ 0.036). There was also a trend toward vaccinated IL6 Ϫ/Ϫ mice having fewer tetramer-positive Blastomyces-specific CD4 ϩ T cells and these cells representing a smaller portion of CD4 ϩ T cells in the lungs of IL6 Ϫ/Ϫ mice than their WT counterparts (Fig. 4C).
Past reports have sought to characterize cellular sources of IL-17 in blastomycosis and have implicated, in addition to Th17 cells, Tc17 cells, ␥␦ T cells, natural Th17 (nTh17) cells, and innate lymphoid cells (18,19). In fact, previous reports have demonstrated that the loss of IL-6 leads to a decrease in Tc17 responses to Blastomyces spp., though this was carried out in CD4 ϩ T cell-deficient mice (19). Though outside the scope of this study, it is possible that other lymphocyte populations in addition to Th17 cells contribute to the IL-17 phenotype we report.
Neutrophils are essential in resistance to B. dermatitidis (20). Compared to vaccinated WT mice, IL6 Ϫ/Ϫ mice evinced a significant paucity of activated neutrophils in the lung, as well as total neutrophils, especially when corrected for fungal burden (Fig. 4D  and S2B). In contrast to the model of acquired resistance, antibody neutralization of IL-6 in unvaccinated WT mice failed to influence resistance to infection (Fig. S2C). Likewise, unvaccinated WT and IL6 Ϫ/Ϫ mice did not differ on any endpoint analyzed (Fig. 4A to  D). Taken together, our findings suggest that IL-6 plays a more prominent role in acquired versus innate resistance to this infection.
Our findings demonstrate that IL-6 is critical to adaptive immunity to B. dermatitidis infection. Considering that IL-6 promotes the development of IL-17-producing CD4 ϩ T cells (17), our findings in both the murine model and the human donors are consistent with past work that identified a role for IL-17 in acquired resistance to experimental blastomycosis in mice (9). While we could not compare memory CD4 ϩ T cell responses to Blastomyces spp. in Hmong and European patients, our use of C. albicans as a surrogate stimulus raises the possibility that lower IL-6 production in the WI Hmong relative to other populations may have broader implications for antifungal immunity.
In summary, we characterize genetic differences that may underlie susceptibility to an endemic mycosis and experimentally analyze population differences in antifungal immunity. Analysis of Hmong blastomycosis patient genomes revealed 113 candidate susceptibility variants. We experimentally validated the IL6 locus, which harbors 25 of Genetic Susceptibility to Blastomycosis ® the candidate variants, more than any other antifungal immunity gene in our analysis. Cells from Hmong donors produced less IL-6 than did cells from European donors and blunted IL-17 production by memory CD4 ϩ T cells in response to fungi. We also formally demonstrated that the loss of IL-6 engenders susceptibility and blunts IL-17producing CD4 ϩ T cells in a murine model of acquired resistance to B. dermatitidis.
Our data in humans and mice support our contention that genetic variation at the IL6 locus, and altered IL-6 responses to fungi, contribute to susceptibility to fungal infection. Loss of IL-6 signaling in humans confers a risk of pneumonia, as evinced by long-term follow-up of patients treated with IL-6 receptor (IL-6R) blockade (21). IL-6 promotes Th17 responses, and a number of genetic defects affecting Th17 have been reported to alter susceptibility to other fungal pathogens, including other endemic mycoses (4) and C. albicans (22). Th17 responses also play a role in the control of other intracellular pathogens, including Mycobacterium tuberculosis (23,24). Intriguingly, the genotype at one candidate susceptibility variant that we report near IL6 is correlated with susceptibility to tuberculosis and hepatitis B virus in Han Chinese populations, with the proposed mechanisms involving the development of Th17 responses (16,25). Clinical recommendations regarding patients on immunologic blockade with biological agents already recommend screening and prophylaxis for granulomatous diseases, including tuberculosis and histoplasmosis, though this is best characterized in the context of tumor necrosis factor alpha (TNF-␣) blockade (26). It may be prudent for clinicians managing these patients to also have higher suspicion of blastomycosis. A number of highly linked candidate variants at the IL6 locus are eQTLs for AS-IL6, a long noncoding RNA that was recently reported to be involved in IL-6 production following stimulation (11). Indeed, one candidate variant is transcribed in AS-IL6. We hypothesize that candidate variant(s) may engender susceptibility to blastomycosis in humans by altering the regulation of IL6, possibly by way of quantitatively or qualitatively altering AS-IL6. Our future experiments will be directed at characterizing the role of select SNPs on the IL-6/AS-IL6 system.
Although we were constrained by the fact that blastomycosis is a relatively uncommon infection, we circumvented this problem by implementing an innovative ROH mapping strategy to investigate the genetic basis of a rare condition and implicate a potential susceptibility locus.
We acknowledge that the size of our study is small and that identification of candidate susceptibility variants was biased toward 108 genes postulated or implicated in immunity to fungi. However, we also created additional criteria for prioritization, specifically, polymorphisms flagged as regulatory, and those that are reported eQTLs. Of 113 candidate variants we report here, only 26 were identified solely on the basis of their proximity to potential antifungal immunity genes. Future studies will seek to apply unbiased informatics approaches, including studies with animal models that are not feasible to carry out with human donors when studying a rare disease. These studies may uncover numerous susceptibility loci not previously implicated in antifungal immunity, which will enable the discovery of polygenic susceptibility variants. We are optimistic that such studies will facilitate future human genetic studies, including broader analysis of genetic susceptibility to blastomycosis among the Hmong people.

MATERIALS AND METHODS
Ethics statement. All human subject research was approved by the institutional review board of the University of Wisconsin-Madison (project ID 2013-1139). Studies involving laboratory animals were approved by the IACUC of UW-Madison (project ID M005891).
Whole-genome sequencing. After obtaining informed consent, we collected saliva samples from individuals of Hmong ancestry. All but one donor had a history of treatment for confirmed blastomycosis, including five individuals who presented with blastomycosis associated with an outbreak in 2013 in Marathon County, WI (1). Illumina sequencing was performed at the University of Wisconsin Biotechnology Center or the Broad Institute. Variant calls for both sets of sequencing runs were done by the Broad Institute using their established best practices pipeline for quality control and mapping of human sequencing data. The GRCh37/hg19 build was used as the human genome reference.
Runs of homozygosity. Prior to identification of ROH, variants were filtered with VCFtools using the following parameters: -max-missing 0.5, -minQ 30, -minDP 3, -remove-indels, -max-alleles 2, -hwe 0.0001. BCFtools/RoH (27) and GARLIC (28) were used to identify ROH in each study participant. Allele frequencies were calculated within GARLIC by resampling 50 individuals from nine unrelated Hmong individuals using the -freq-only flag. The parameters used to calculate ROH with Garlic were garlic -build hg19 -error 0.001 -auto-winsize -auto-overlap-frac -winsize 50; the -freq-file generated from the unrelated individuals was provided; and default clusters as well as a run with -nclust 5 were run. The same frequencies were provided to the BCFtools/RoH algorithm: bcftools roh -G30 -AF-file. Visualization of ROH identified by both programs was performed with ggplot2 in R (https://cran.r-project.org/web/ packages/ggplot2/index.html). Recognizing imprecision in identifying ROH, we merged individual ROH calls from the two algorithms using bedtools merge. We looked for overlapping ROH among cases using BEDOPS (29) and focused on regions in which at least seven of the nine cases had an ROH identified; we refer to these regions as "consensus ROH." Candidate susceptibility variants. Variants within consensus ROH that were fixed in Hmong blastomycosis cases and found at Ͻ20% frequency in the European (EUR) Superpopulation of the 1000 Genomes Project (30) were inspected to identify candidate susceptibility variants. We focused on variants that fulfilled at least one of three criteria, (i) proximity to antifungal immunity genes, (ii) annotation as regulatory, or (iii) reported expression quantitative trait loci (eQTL).
Starting from a list of genes in an antifungal immunity array and adding genes based on the literature, we interrogated 108 genes known or suspected to be involved in human immune responses to fungi (Table S2). We flagged SNPs occurring within 100 kb on either side of the transcriptional start site of these genes. SNPs within consensus ROH were annotated using the Ensembl Variant Effector Predictor (31) and cross-referenced with the 1000 Genomes and Genotype-Tissue Expression (GTEx; Statistical analysis. Data were analyzed in GraphPad Prism 7.04 using the two-tailed Mann-Whitney test. A P value of Յ0.05 was defined as statistically significant.