Systems genetics uncover new loci containing functional gene candidates in Mycobacterium tuberculosis-infected Diversity Outbred mice

Abstract

Mycobacterium tuberculosis infects two billion people across the globe, and results in 8-9 million new tuberculosis (TB) cases and 1-1.5 million deaths each year. Most patients have no known genetic basis that predisposes them to disease. Here, we investigate the complex genetic basis of pulmonary TB by modelling human genetic diversity with the Diversity Outbred mouse population. When infected with M. tuberculosis, one-third develop early onset, rapidly progressive, necrotizing granulomas and succumb within 60 days. The remaining develop non-necrotizing granulomas and survive longer than 60 days. Genetic mapping using immune and inflammatory mediators; and clinical, microbiological, and granuloma correlates of disease identified five new loci on mouse chromosomes 1, 2, 4, 16; and three known loci on chromosomes 3 and 17. Further, multiple positively correlated traits shared loci on chromosomes 1, 16, and 17 and had similar patterns of allele effects, suggesting these loci contain critical genetic regulators of inflammatory responses to M. tuberculosis. To narrow the list of candidate genes, we used a machine learning strategy that integrated gene expression signatures from lungs of M. tuberculosis-infected Diversity Outbred mice with gene interaction networks to generate scores representing functional relationships. The scores were used to rank candidates for each mapped trait, resulting in 11 candidate genes: Ncf2, Fam20b, S100a8, S100a9, Itgb5, Fstl1, Zbtb20, Ddr1, Ier3, Vegfa, and Zfp318. Although all candidates have roles in infection, inflammation, cell migration, extracellular matrix remodeling, or intracellular signaling, and all contain single nucleotide polymorphisms (SNPs), SNPs in only four genes (S100a8, Itgb5, Fstl1, Zfp318) are predicted to have deleterious effects on protein functions. We performed methodological and candidate validations to (i) assess biological relevance of predicted allele effects by showing that Diversity Outbred mice carrying PWK/PhJ alleles at the H-2 locus on chromosome 17 QTL have shorter survival; (ii) confirm accuracy of predicted allele effects by quantifying S100A8 protein in inbred founder strains; and (iii) infection of C57BL/6 mice deficient for the S100a8 gene. Overall, this body of work demonstrates that systems genetics using Diversity Outbred mice can identify new (and known) QTLs and functionally relevant gene candidates that may be major regulators of complex host-pathogens interactions contributing to granuloma necrosis and acute inflammation in pulmonary TB.

Fig 1. Mouse survival following exposure to a low dose of aerosolized M. tuberculosis strain Erdman.

Diversity Outbred (DO) mice (n = 680, brown solid line), and the eight inbred founder strains (n = 15 to 78, colored lines) were infected with M. tuberculosis strain Erdman bacilli by aerosol. Panel A: Cumulative survival. Approximately 30% of the DO population succumbed to pulmonary TB by 60 days post infection and approximately 70% succumbed to between 60- and 600-days post infection. Of the eight inbred founder strains, survival studies were completed for the C57BL/6J inbred strain; the other seven inbred founder strains were euthanized 40 days after M. tuberculosis infection. No inbred founder or non-infected (NI) DO mice (n = 53, dashed line) showed mortality in that period. Panel B: Daily mortality of M. tuberculosis-infected DO mice, highlighting the early wave of mortality that peaked between 25- and 35-days post infection.

Fig 2. Diversity Outbred (DO) mice develop a spectrum of histopathological lung lesions due to M. tuberculosis infection.

Lung lobes were formalin-fixed, paraffin-embedded, sectioned, and stained with hematoxylin & eosin. Panel A: Lung section from a non-infected DO mouse. Panels B through I: Lung sections from M. tuberculosis-infected DO mice euthanized eight weeks after infection show a spectrum of lesions from mild to severe (upper left to bottom right); focal lesions (e.g., Panel B) to diffuse infiltration (Panel H); and include necrotizing (Panels C, E, F, I) and non-necrotizing granulomas (Panels B, D, G, H). Low magnification (15X).

Fig 3. M. tuberculosis infection of Diversity Outbred (DO) mice induces clinical, microbiological, and immune/inflammatory phenotypes that are positively and significantly correlated.

(A) Boxplots of phenotypes in non-infected (n = 77) and infected (n = 930) mice at euthanasia, which ranged from 5 to 553 days post-infection. Each box shows the median (central line) and the inter-quartile range in the box. The whiskers show 90% of the data and the dots beyond the whiskers show data points outside of the 90th percentile of the distribution. Significant differences at a Bonferroni corrected p-value of 0.004 are shown with an asterisk. (B) Pearson correlation of selected phenotypes with n = 930 infected mice. Each cell in the upper triangle contains an ellipse in which the color and shape show the strength of the correlation. Positive correlations are shown in blue and negative correlations are shown in red. Strong correlations are shown by narrow ellipses, and weak correlations are shown by round ellipses. The more circular the ellipse, the closer the correlation is to zero. Numerical values of the Pearson correlation are shown in the lower triangle and range from -0.14 to 0.84. Phenotypes are ordered by hierarchical clustering of the similarity in correlations. Absolute Pearson correlations greater than 0.105 are significant at a Bonferroni-corrected p-value of 0.0007 and are shown with an asterisk.

Fig 4. The heatmap of linkage mapping peaks shows patterns of common genetic regulation.

Quantitative trait locus (QTL) mapping results from M. tuberculosis infected Diversity Outbred mice that were euthanized on or prior to 250 days post infection (if morbidity developed), and complete genotype and phenotypic trait data which also passed QA/QC (n = 853). The mouse genome, from chromosome 1 through X, is shown on the horizontal axis. Phenotypes are shown on the vertical axis. Each cell shows the logarithmic of the odds (LOD) score on one chromosome for the phenotype listed on the left, colored by the color scale. Higher LOD scores (oranges and yellows) meet or exceed statistical thresholds and show the significant genotype-phenotype associations. The phenotypes are hierarchically clustered based on the correlation between LOD curves, i.e., phenotypes with similar LOD curves are clustered next to each other.

Fig 5. Overview of gene prioritization methods.

Traits were mapped to identify positional candidate genes in QTLs. Gene expression data were analyzed for differential gene expression. The top 500 differentially expressed genes (DEG) were used to train SVMs to distinguish these trait-related genes from other genes in the genome using the FNTM mouse lung network. The fitted models were used to score positional candidates in each trait QTL. Positional candidates were then ranked as trait-related based on their functional scores.

Fig 6. Quantitative Trait Locus (QTL) mapping results of first principal component (PC1) of CXCL1, CXCL2, M. tuberculosis burden, and MMP8 identifies Dots1 on chromosome 1, containing the gene candidates Fam20b and Ncf2.

Panel A: LOD curve for PC1 between 150 and 160 Mb on chromosome 1 with peak near 155.36 Mb. Panel B: Founder allele effects for PC1 in the same genomic interval. Each colored line is the best linear unbiased predictor for one of the founder alleles. Founder colors are shown in the upper left. Panel C: LOD score of the imputed SNPs in the same genomic interval. Each point represents the LOD score of one imputed SNP. Panel D: Genes in the confidence interval. Panel E: Functional scores for genes in the chromosome 1 QTL. Each dot represents a single gene. Its position on the x axis is its position within the QTL. Its position on the y axis is the functional -log10(UPPR) derived from the SVM. Points are colored based on is correspondence with the trait—green with CXCL1, orange with MMP8, and blue with M. tuberculosis burden. Fam20b and Ncf2 genes had the highest functional scores. QTL mapping results from n = 853 M. tuberculosis infected Diversity Outbred mice up to 250 days post infection.

Fig 7. Quantitative Trait Locus (QTL) mapping of lung S100A8 identifies Dots3 on chromosome 3, containing the gene candidates S100a8 and S100a9.

Panel A: LOD score in the confidence interval from 90 to 92 Mb on chromosome 3. Panel B: Founder allele effects within the confidence interval. Panel C: LOD score of the imputed SNPs with the highest SNPs colored in red. Panel D: Genes in the same interval. The gene S100a8 is directly under the SNPs with the highest LOD scores. Panel E: Functional scores for genes in the chromosome 3 QTL. Each dot represents a single gene. Its position on the x axis is its position within the QTL. Its position on the y axis is the functional -log10(UPPR) derived from the SVM. The gene S100a8 had the highest functional score. QTL mapping results from n = 853 M. tuberculosis infected Diversity Outbred mice up to 250 days post infection.

Fig 8. Quantitative Trait Locus (QTL) mapping results of first principal component (PC1) of M. tuberculosis burden, weight loss, and granuloma necrosis identifies Dots5 on chromosome 16, containing gene candidates Fstl1, Itgb5, and Zbtb20.

M. tuberculosis lung burden, weight loss, and granuloma necrosis map to a region on chromosome 16 near 38 Mb. Panel A: LOD profile for M. tuberculosis in the confidence interval. The genomic position on chromosome 16 is on the horizontal axis and the LOD score is on the vertical axis. Panel B: Founder allele effect in the confidence interval. The vertical axis shows the estimates effect of gaining one founder allele. Panel C shows the SNP LOD score for association mapping using imputed SNPs. Panel D: Genes in the confidence interval. Panel E: Functional scores for genes in the chromosome 16 QTL. Each dot represents a single gene. Its position on the x axis is its position within the QTL. Its position on the y axis is the functional -log10(UPPR) derived from the SVM. Points are colored based corresponds to the trait—green with weight loss and blue with M. tuberculosis burden. Fstl1 was the top ranked gene overall followed by Itgb5 and Zbtb20. QTL mapping results from n = 853 M. tuberculosis infected Diversity Outbred mice up to 250 days post infection.

Fig 9. Quantitative Trait Locus (QTL) mapping results of first principal component (PC1) of granuloma necrosis, M. tuberculosis burden, weight loss, CXCL1 and MMP8 identifies Dots7 on chromosome 17, which contains gene candidates Vegfa, Ier3, Ddr1, and Zpf318.

Panel A: LOD score for PC1 of lung granuloma necrosis ratio, M. tuberculosis burden, MMP8, CXCL1, and weight loss in the interval where the phenotypes map. Panel B: Founder allele effects for the two peaks. Panel C: LOD score of the imputed SNPs in the interval, with the highest scoring SNPs colored in red. Panel D: Functional scores for genes in the chromosome 17 QTL. Each dot represents a single gene. Its position on the x axis is its position within the QTL. Its position on the y axis is the functional -log10(UPPR) derived from the SVM. Points are colored based on correspondence to the traits—blue with CXCL1 and green with granuloma necrosis. Vegfa was the top ranked gene overall followed by Ier3, Ddr1, and Zpf318. QTL mapping results from n = 853 M. tuberculosis infected Diversity Outbred mice up to 250 days post infection.