Stratified analyses refine association between TLR7 rare variants and severe COVID-19

Abstract

Despite extensive global research into genetic predisposition for severe COVID-19, knowledge on the role of rare host genetic variants and their relation to other risk factors remains limited. Here, 52 genes with prior etiological evidence were sequenced in 1,772 severe COVID-19 cases and 5,347 population-based controls from Spain/Italy. Rare deleterious TLR7 variants were present in 2.4% of young (<60 years) cases with no reported clinical risk factors (n = 378), compared to 0.24% of controls (odds ratio [OR] = 12.3, p = 1.27 × 10^-10). Incorporation of the results of either functional assays or protein modeling led to a pronounced increase in effect size (OR_max = 46.5, p = 1.74 × 10^-15). Association signals for the X-chromosomal gene TLR7 were also detected in the female-only subgroup, suggesting the existence of additional mechanisms beyond X-linked recessive inheritance in males. Additionally, supporting evidence was generated for a contribution to severe COVID-19 of the previously implicated genes IFNAR2, IFIH1, and TBK1. Our results refine the genetic contribution of rare TLR7 variants to severe COVID-19 and strengthen evidence for the etiological relevance of genes in the interferon signaling pathway.

Keywords: SARS-CoV-2; burden analysis; host genetics; immune deficiency; infection; innate immunity; rare variants; targeted sequencing; toll-like receptor 7; variant collapsing analysis.

Figure 1

Study design and cohort characteristics (A) Candidate genes included in targeted sequencing, grouped according to source of evidence (details in Table S1). Genes known to cause human inborn errors of immunity are highlighted in bold, and genes excluded during quality control due to low sequencing coverage are crossed out. (B) Workflow describing the main steps of sample preparation, genotyping, sequencing, and computational processing. Boxes colored in gray indicate steps that were performed in previous studies.^, MIP, molecular inversion probe; PCA, principal-component analysis; QC, quality control. (C) Number of individuals in the Italian (left) and Spanish (right) subcohorts. The number of COV_hosp individuals with no reported risk factors (as described in Table 1) is highlighted in red. The proportion of females is shown in parentheses. (D) Number of variants observed in the cohort in relation to their minor allele frequency (MAF). In the present study, variants with MAF <1% were denoted as rare variants, while all others were considered common. Intensity of color shading indicates whether (dark) or not (light) variants have been reported with allele frequency in gnomAD r2.1 exomes.

Figure 2

Association analysis for individual variants The p values (y axis, negative log10) obtained in the association analysis of 1,211 non-singleton variants from the POP_all analysis. Variants are grouped according to the genes (x axis, sorted alphabetically) in which they are located. Results for case-control definitions other than POP_all are provided in Figure S3. Dotted line: Lenient significance threshold, correcting for the number of variants tested (α = 4.1 × 10⁻⁵). Dashed line: Strict significance threshold, also taking into account multiple testing due to additional case-control definitions (α = 6.7 × 10⁻⁶). Variants with p values below the lenient significance threshold are marked in green and were only found in genes selected based on prior GWAS evidence, i.e., FYCO1 and XCR1 at 3p21.31, IFNAR2 at 21q22.11.

Figure 3

Results of the gene-based collapsing analysis for rare variants with MAF <0.1% The p values (y axis, negative log10) are plotted for 52 genes (x axis, sorted alphabetically). The various case-control definitions (see Table 1; excluding sex-stratified analyses) are depicted as symbols, while variant deleteriousness classes are coded according to color (M1: pLoF, M3 and M4: pLoF and moderate effect variants including missense in two graduations, C10+M1: CADD >10 or pLoF, SYN: synonymous, see subjects and methods). Dashed line: Strict significance threshold, correcting for all tests conducted: (α = 8.7 × 10⁻⁶). Dotted line: Lenient significance threshold, correcting for the number of genes tested (α = 9.6 × 10⁻⁴). Results for sex-stratified analyses and variants with MAF <1% are provided in Figure S4.

Figure 4

Forest plot for TLR7 rare variant gene burden according to variant classification (A) Odds ratios (ORs) of collapsed variants in TLR7 are shown for POP_lowrisk at different minor allele frequency groups (MAF) and deleteriousness predictions (class). Within each group, results are presented for all individuals and for sex-stratified analyses. Error bars indicate 95% confidence intervals. SYN, synonymous; C10+M1, CADD>10 or pLoF; bLoF, biochemical evidence for a loss-of-function effect; 3D-P, variant class based on 3D protein structure, see subjects and methods. SYN variants with MAF <0.001 were only present in controls (OR = 0.0, no confidence interval calculable). (B) Presence of 3D-P TLR7 (MAF<0.1%) variant carriers (black dots) in all COV_hosp individuals (gray blue), COV_hosp with age <60 y (light blue) and COV_hosp with no reported risk factors (“no RF,” dark blue). The number of individuals within each set is indicated by area and is specified in the outer legend. Percentages in brackets represent carrier ratios.

Figure 5

Location of rare TLR7 variants within TLR7 protein domains (A) Rare, deleterious TLR7 variants (POP_lowrisk, C10+M1, MAF <0.1%) are mapped on the protein domains of TLR7 (x axis: amino acid position). Phenotype, according to the POP_lowrisk case-control definition, and the sex of variant carriers is indicated by color or caption. Variants of carriers previously reported in Asano et al. (see subjects and methods and Table S2) are indicated by asterisks (∗). TLR7 domains: LRR-NT (leucine-rich repeat, N terminal, aa 27–65) orange; LRR regions 1–26 (aa 66–786) blue; LRR-CT (leucine-rich repeat, C terminal, aa 787–839), violet; TIR (Toll/interleukin-1 receptor) domain (aa 889–1033), green. (B) TLR7 dimer overview, interface highlighted as gray surface (also in C). (C) Non-synonymous variants from (A) are highlighted in the 3D conformation of one TLR7 subunit (PDB ID: 5GMH) and are presented from two angles. Phenotype (POP_lowrisk, see A) and sex of the variant carriers are indicated by color coding. Variants within 5 Ångström of the subunit interface are highlighted by a hash (#, also in A). Variants located downstream of position T858 could not be plotted due to absence of the respective residues from the structure. Visualized using PyMOL Molecular Graphics System (Version 2.5.5 Schrödinger, LLC).