Genetic Analysis of Bphse: a Novel Gene Complementing Resistance to Bordetella pertussis-Induced Histamine Sensitization

Histamine is a bioactive amine associated with a plethora of normal and pathophysiological processes, with the latter being dependent on both genetic and environmental factors including infectious agents. Previously, we showed in mice that susceptibility to Bordetella pertussis and pertussis toxin (PTX) induced histamine sensitization (Bphs) is controlled by histamine receptor H1 (Hrh1/HRH1) alleles. Bphs susceptible and resistant alleles (Bphss/Bphsr) encode for two-conserved protein haplotypes. Given the importance of HRH1 signaling in health and disease, we sequenced Hrh1 across an extended panel of laboratory and wild-derived inbred strains and phenotyped them for Bphs. Unexpectedly, eight strains homozygous for the Bphsr allele phenotyped as Bphss, suggesting the existence of a modifying locus segregating among the strains capable of complementing Bphsr. Genetic analyses mapped this modifier locus to mouse chromosome 6; designated Bphs-enhancer (Bphse), within a functional linkage disequilibrium domain encoding multiple loci controlling responsiveness to histamine (Bphs/Hrh1 and Histh1-4). Interval-specific single-nucleotide polymorphism (SNP) based association testing across 50 laboratory and wild-derived inbred mouse strains and functional prioritization analyses resulted in the identification of candidate genes for Bphse within a ∼5.5 Mb interval (Chr6:111.0-116.4 Mb), including Atg7, Plxnd1, Tmcc1, Mkrn2, Il17re, Pparg, Lhfpl4, Vgll4, Rho and Syn2. Taken together, these results demonstrate the power of combining network-based computational methods with the evolutionarily significant diversity of wild-derived inbred mice to identify novel genetic mechanisms controlling susceptibility and resistance to histamine shock.


Histamine (2-[4-imidazole]-ethylamine; HA) is an endogenous biogenic
monoamine that is synthesized, stored intracellularly within granules; following cellular activation HA is released by mast cells, basophils, platelets, neurons, and enterochromaffin-like cells in the stomach [1]. After release, free HA mediates its pleiotropic effects by binding to four different seven-transmembrane G-protein-coupled receptors (GPCRs): histamine receptor H1-H4 (HRH1, HRH2, HRH3, and HRH4), differentially expressed on target cells in various tissues [2]. HA acting through these receptors influences a diverse array of physiological processes, including brain function, neurotransmission, secretion of pituitary hormones, cell proliferation and differentiation, hematopoiesis, embryonic development, wound healing and regeneration, and the regulation of gastrointestinal, cardiovascular, and secretory functions [2]. In addition, HA plays a major role in inflammation and the regulation of innate and adaptive immune responses in both normal and pathologic states [3,4].
Historically, HA is most well-known for its role in shock and anaphylaxis [5]. It was first isolated from the parasitic mold ergot of rye (Claviceps purpurea) and then synthesized by the decarboxylation of histidine [6][7][8]. HA was shown to elicit anaphylactic shock-like symptoms when injected into mammals, including bronchiolar constriction, constricted cardiac and pulmonary arteries, and stimulated cardiac contraction [9][10][11].
Further research firmly established HA as a natural constituent of the body and a mediator of anaphylactic shock [5]. There is significant variability in susceptibility to HA-shock among animal species, with guinea pigs and rabbits being highly susceptible, versus mice and rats which are generally remarkably tolerant to in vivo injections of HA [12].
Interestingly, prior exposure to Bordetella pertussis (B. pertussis) or purified pertussis toxin (PTX) overcomes HA resistance among a subset of laboratory derived inbred strains of mice [13]. This phenotype is designated Bphs for B. pertussis-induced HA sensitization [14,15]. Bphs-susceptible (Bphs s ) strains die within 30 minutes following HA injection; which is thought to result from hypotensive and hypovolemic shock while Bphs-resistant (Bphs r ) strains remain healthy [16]. The sensitizing activity elicited by exposure to B.
Our previous genetic studies mapped the autosomal dominant Bphs locus controlling susceptibility to mouse chromosome 6 (Chr6) and identified it as histamine receptor H1 (Hrh1/HRH1) [15,20]. Susceptibility segregates with two conserved Hrh1/HRH1 haplotypes, mice with the Hrh1 s /HRH1 s allele (encoding Pro 263 , Val 312 , Pro 330 ) phenotype as Bphs s while mice with the Hrh1 r /HRH1 r allele (encoding Leu 263 , Met 312 , Ser 330 ) are Bphs r . These amino acid changes occur within the third intracellular loop of this G protein-coupled receptor (GPCR): a domain implicated in signal transduction, protein folding, and trafficking. Functionally, HRH1 s and HRH1 r alleles equally activate Gαq/11, the G protein family members that couple HRH1 signaling to second messenger signaling pathways, indicating that the genetic control of susceptibility and resistance to Bphs is not inherently due to differential activation of either Gαq or Gα11 [21]. However, the two alleles exhibit differential cell surface expression and altered intracellular trafficking, with the HRH1 r allele selectively retained within the endoplasmic reticulum

RESULTS
Hrh1/HRH1 alleles are highly conserved in mice. We undertook a genetic approach to screen for evolutionarily selected mechanisms that may be capable of modifying the Bphs r phenotype. Toward this end, we sequenced ~500 bp stretch of genomic DNA encompassing the third intracellular loop of Hrh1/HRH1 across 91 laboratory and wild-derived inbred strains of mice (Table 1). Surprisingly, other than the three amino acid changes (Pro 263 , Val 312 , Pro 330 → Leu 263 , Met 312 , Ser 330 ) described earlier as "susceptible" (Hrh1 s /HRH1 s ) and "resistant" (Hrh1 r /HRH1 r ) haplotypes [15], we did not identify any additional non-synonymous amino acid changes. Of the 91 strains, 22 carry the Hrh1 r /HRH1 r allele, whereas 69 carry the Hrh1 s /HRH1 s allele ( Table 1). We next mapped the evolutionary distribution of the two alleles onto a mouse phylogenetic tree (Supplementary Figure 1) [22]. The Hrh1 r /HRH1 r allele was primarily restricted to wildderived group 7 strains and a select sub-branch of group 1 Bagg albino derivatives, whereas the Hrh1 s /HRH1 s allele was distributed across groups 2-6, which encompasses Swiss mice, Japanese and New Zealand inbred strains, C57/58 strains, Castle mice, C.C.
Little DBA and related strains.  [28][29][30]. This genetic variability represents a rich source of evolutionarily selected diversity and has the potential to lead to the identification of genes controlling novel regulatory features arising from host-pathogen co-evolutionary adaptations.
To screen for functional modifying loci capable of complementing Hrh1 r /HRH1 r , we phenotyped a panel of group 1 (Bagg albino derivatives) and 7 (wild-derived) mice that genotyped as Hrh1 r /HRH1 r for susceptibility to Bphs. While nine Hrh1 r /HRH1 r strains tested were Bphs r as expected, we found eight that were remarkably susceptible to Bphs ( Table 2). Importantly, these Bphs s strains are confined primarily to group 7 wild-derived strains (Supplementary Figure 1) in contrast to Hrh1 r /HRH1 r strains from Group 1 that were mostly Bphs r . Moreover, comparison of the entire Hrh1 gene (Chr6:114,397,936-114,483,296 bp) between several of the Group 1 and Group 7 phenotyped strains found no segregating non-synonymous structural variants (data not shown) between Bphs s and Bphs r strains suggesting that the complimenting locus/loci in Group 7 is independent of additional previously unidentified Hrh1 s /HRH1 s structural variants.
To confirm the existence of a modifying locus capable of restoring Bphs s in mice with a Hrh1 r /HRH1 r allele and to assess its heritability, we selected a subset of Bphs s -Hrh1 r /HRH1 r (MOLF, PWK) and Bphs r -Hrh1 r /HRH1 r (AKR, CBA, C3H, MRL) strains for follow-up studies. We studied F1 hybrids between the select strains of interest and HRH1knockout B6 mice (HRH1KO), which lack a functional Hrh1 s /HRH1 s gene required for   (Figure 2A and Supplementary Table 2).
There was no difference in predicted candidate genes using either the smaller dataset (13,257 SNPs)  Each of the seven gene sets define a putative Bphs-related process that forms a distinct subnetwork of the full functional genomic network. Using this approach, we identified several hundred genes within the Bphse congenic locus that are functionally associated with each biological process, and thus could be gene candidates (Supplementary Table 4).
Genes that are predicted to be highly functionally related to a trait may not have functional variant alleles segregating in the study population, and may therefore be unlikely to drive the observed strain differences in Bphs s . Using the list of polymorphic genes identified through high-resolution genetic association testing (Figure 2), we normalized and plotted the respective genetic association score (-log10 pEMMA) with functional enrichment (-log10 FPR) to focus on genes that overlap both approaches ( Figure 3A). The final ranking was calculated by defining a final gene score (Scg) for each gene, which is the sum of the (normalized) -log10(FPR) and the -log10 (pEMMA) ( Figure 3B).  Figure 1) and replicated the phenotype. To our knowledge, this is the first study assessing Bphs s in multiple wild-derived inbred strains of mice, and clearly establish their utility in identifying novel genetic mechanisms controlling HA-shock.
Aside from genetics, several factors could influence HA sensitivity after exposure to B. pertussis and PTX including age, sex, and route of sensitization/challenge [12]. In our phenotyping experiments, we used 8-12-week-old mice of each sex and did not find any sex differences. This agrees with earlier studies that found no sex differences in Bphs s [47]. We also tested the route of administration of PTX and HA challenge using the intraperitoneal and intravenous routes and found no difference (data not shown). We have not tested the effect of age on Bphs s amongst the various strains; however, work from Munoz and others have reported a significant effect of age [12]. It is possible that some of the strains that are Bphs r will exhibit Histh s as they age or following treatment with  (Figure 2).
Recently, a quantitative trait gene prediction tool has been described that utilizes functional genomics information (gene co-expression, protein-protein binding data, ontology annotation and other functional data) to rank candidate genes within large QTLs associated with a respective phenotype [43]. This methodology uses biological prior knowledge to predict candidate genes that could influence multiple pathways affecting the phenotype. We utilized this approach for Bphs, which is known to involve cardiac, vascular, and anaphylactic mechanisms [44,45]. Because the selection of phenotypeassociated gene sets is critical for final gene predictions, several terms were used to incorporate sub-phenotypes equivalent to Bphs in the expectation that use of multiple terms would help identify candidate loci for Bphse. Integration of functional predictions with genetic association (Scg, Figure 3) allowed us to focus on only those candidates that reached significance in both approaches.
Given that there is differential cell surface expression of HRH1 depending on the haplotype, it is tempting to speculate that Bphse may aid in the folding, trafficking and/or Importantly, the fact that Bphse resides within a smaller functional LD that includes DNA isolation and genotyping. DNA was isolated from mouse tail clippings as previously described [13]. Briefly, individual tail clippings were incubated with cell lysis buffer (125 mg/ml proteinase K, 100 mM NaCl, 10 mM Tris-HCl (pH 8.3), 10 mM EDTA, 100 mM KCl, 0.50% SDS, 300 ml) overnight at 55 o C. The next day, 6M NaCl (150 ml) was added followed by centrifugation for 10 min at 4 o C. The supernatant layer was transferred to a fresh tube containing 300 µl isopropanol. After centrifuging for 2 min, the supernatant was discarded, and the pellet washed with 70% ethanol. After a final 2 min centrifugation, the supernatant was discarded, and DNA was air dried and resuspended in TE. Genotyping was performed using microsatellite, sequence specific, and Hrh1 primers (Supplemental Table 5).

Low-resolution interval-specific targeted genetic association testing.
Genotype data (SNPs in both coding and non-coding) of 50 mouse strains Trait-related gene sets. The positional candidate genes were ranked based on their predicted association with seven functional terms related to the Bphs phenotype: "Cardiac", "G-protein coupled receptor", "Histamine", "Pertussis toxin", "Type I hypersensitivity", "Vascular Permeability", and "ER/EMC/ERAD". Gene Weaver [72] was used to identify genes annotated with each term. Each term was entered the Gene Weaver homepage (https://geneweaver.org)and search restricted to human, rat, and We used 10-fold cross validation and a linear kernel. We also trained each SVM over a series of cost parameters identified by iteratively narrowing the cost parameter window to identify a series of eight cost parameters that maximized classification accuracy. We then used the training modules to score each positional candidate gene in the Bphse locus.
To compare scores across multiple trained models, we converted SVM scores to false positive rates.
Combined gene score. To create the final ranked list of positional candidate genes, we combined the SNP association scores with the functional predictions derived from the SVMs. We scaled each of these scores by its maximum value across all positional candidates and summed them together to derive a combined gene score (Scg) that incorporated both functional predictions and genetic influence: where the denominators of the two terms on the right-hand side are the maximum values of -log10(pEMMA) and -log10(FPRSVM) over all positional candidates in Bphse, respectively, which normalizes the functional and positional scores to be comparable to each other. SNPs were assigned to the nearest gene within 1Mb. If more than one SNP was assigned to a gene, we used the maximum negative log10 p value among all SNPs assigned to the gene.          This list was generated using a stringent cut-off (p<3.81E-06) and a moderate cut-off (p<5.00E-02) using Efficient Mixed Model Association (EMMA). Only the most significant SNP tagging a gene is shown. The functional location of variants is listed using the following notation: intronic (I), synonymous change (Cs), 3' untranslated region variant (U3), 5' untranslated region variant (U5), not characterized (NC) and no gene (NA). 5.45E-05 This list was generated using a moderate cut-off (p<5.00E-02) using Efficient Mixed Model Association (EMMA). Only the most significant SNP tagging a gene is listed. The functional location of variants is listed using the following notation: intronic (I), synonymous change (Cs), nonsynonymous change (Cn), 3' untranslated region variant (U3), 5' untranslated region variant (U5), not characterized (NC) and no gene (NA).

Supplementary Table 4.
List of genes predicted to be functionally associated with Bphs physiological processes ranked by negative log of false positive rate.