A Genome-Wide Association Study of Respiratory Syncytial Virus Infection Severity in Infants

Abstract

Background

Respiratory syncytial virus (RSV) is a significant cause of infant morbidity and mortality worldwide. Most children experience at least one 1 RSV infection by the age of two 2 years, but not all develop severe disease. However, the understanding of genetic risk factors for severe RSV is incomplete. Consequently, we conducted a genome-wide association study of RSV severity.

Methods

Disease severity was assessed by the ReSVinet scale, in a cohort of 251 infants aged 1 week to 1 year. Genotyping data were collected from multiple European study sites as part of the RESCEU Consortium. Linear regression models were used to assess the impact of genotype on RSV severity and gene expression as measured by microarray.

Results

While no SNPs reached the genome-wide statistical significance threshold (P < 5 × 10⁻⁸), we identified 816 candidate SNPs with a P-value of <1 × 10⁻⁴. Functional annotation of candidate SNPs highlighted genes relevant to neutrophil trafficking and cytoskeletal functions, including LSP1 and RAB27A. Moreover, SNPs within the RAB27A locus significantly altered gene expression (false discovery rate, FDR P < .05).

Conclusions

These findings may provide insights into genetic mechanisms driving severe RSV infection, offering biologically relevant information for future investigations.

Respiratory syncytial virus (RSV) causes seasonal epidemics of respiratory infections and is the most common cause of respiratory infection in neonates and young children worldwide. While most experience mild symptoms, for approximately 1 in 56 healthy term infants, infection can lead to hospitalization, placing a significant burden on health care systems, particularly in low- and middle-income settings [1]. In severe cases, the infection can progress to the lower respiratory tract, causing bronchiolitis, pneumonia, and potentially fatal respiratory failure [2]. Identifying genetic risk factors for severe respiratory distress is therefore crucial for managing disease.

RSV is an enveloped, single-stranded, negative-sense RNA virus belonging to the subfamily Pneumoviridae. Its genome encodes 11 structural and nonstructural proteins, essential for viral replication and immune evasion [3]. Infants and children with underlying medical conditions, such as preterm birth, chronic lung disease, or congenital heart disease, are particularly susceptible to severe RSV infections [4, 5]. However, even without underlying risk factors, RSV infection can cause serious disease, and most hospital admissions are among previously healthy infants [6]. Current prevention strategies focus on passive immunization with monoclonal antibodies to high-risk infants during the RSV season or maternal vaccination [7, 8].

Early twin studies suggested a genetic component to developing severe RSV disease, as monozygotic twins demonstrated a higher concordance in hospitalization rates as compared with dizygotic twins [9]. A candidate gene association study (CGAS) conducted in 2007 also identified several innate immune genes (IFNA13, IL15, TLR8, and STAT1) and chemotaxis genes (CCL8, ITGB2, and VCAM1) associated with RSV and influenza [10]. A recent genome-wide association study (GWAS) of RSV bronchiolitis and infection susceptibility highlighted novel candidate single-nucleotide polymorphisms (SNPs), including those within the CXCR3 gene locus, and genes within olfactory and taste receptors [11, 12]. Despite these findings, there is limited functional characterization of genetic risk loci.

The majority of GWAS variants reside in noncoding regions of DNA and are thought to exert their effects through gene transcription regulation. Expression quantitative trait locus (eQTL) analysis can model the effect of genotype on gene expression, highlighting potential causal variants [13]. No eQTL analysis of RSV infection severity has been conducted to date. Consequently, this study has 2 aims: first, to identify genes and pathways associated with the severity of RSV infection through a GWAS; second, to perform eQTL analysis of SNPs associated with RSV severity.

METHODS

Study Design

For our GWAS of RSV severity, we included participants who were recruited as part of the RESCEU project (Respiratory Syncytial Virus Consortium in Europe), which has been described in detail (ClinicalTrials.gov, NCT03756766) [14]. The study prospectively enrolled RSV-infected infants aged <1 year who were otherwise healthy (N = 500) from hospitals and community clinics across the Netherlands, Spain, and United Kingdom during the 2017–2020 RSV seasons (from October to April).

Infants with preexisting medical conditions, including preterm birth (<37 weeks’ gestation), were included as they were suspected to develop severe RSV (n = 50). Exclusion criteria were previous receipt of antivirals, steroids, or montelukast within 7 days preenrollment and/or exposure to an investigational RSV vaccine or medication. From the recruited cohort, 261 cases were included for further quality control and downstream analysis of RSV severity. Samples were selected by clinical data completion and sample availability (Table 1).

Infants were tested for RSV either through point-of-care testing with the Alere i RSV assay (Abbott) in community and hospital settings or through routine tests such as rapid antigen detection or polymerase chain reaction at a central laboratory, where available [15]. RSV disease severity was quantified at the time of recruitment through an extensive clinical questionnaire meant to calculate a ReSVinet score. The score quantifies RSV disease severity on a range of 0 to 20, accounting for respiratory rate, blood oxygen, apnea, dyspnea, feeding intolerance, and other markers [16]. Additional data consisted of hospitalization status and respiratory support requirements (Table 1).

After presentation with RSV disease in the hospital or at health centers, consent for genetic screening was acquired from parents according to the principles of the International Council for Harmonisation’s good clinical practice and national clinical trial regulations. UK study authorization was granted from the Health Research Authority (Integrated Research Application System project 231136), and Research Ethics Committee approval was granted by South Central and Hampshire A (National Research Ethics Service 17/SC/0522). In Spain, this was granted by the Comité de Ética de la Investigación de Santiago-Lugo (2017/395) and in the Netherlands by the Medical Ethical Committee, University Medical Center Utrecht (17/563).

Genome-Wide Association Study

DNA Extraction

DNA was extracted from a buccal swab taken at the time of recruitment and stored in a DNA preservation buffer. Buccal swabs were stored at 20 °C until required for extraction. Buccal swabs were compressed into the medium to ensure suspension of cheek cells. DNA was extracted from buccal swab medium by the QIAGEN DNeasy blood and tissue kit (No. 69504) with a QIASymphony machine to perform the extraction steps according to the standard DNeasy tissue extraction protocol.

SNP Genotyping Assay

All samples were genotyped on the Axiom PMRA platform from Affymetrix with a protocol recommended by the manufacturer.

SNP Quality Control

Prior to association analysis, SNPs with high levels of missingness (>10%) were removed. Likewise, samples with a high proportion of missing data across all genotyped SNPs were removed (>10%), followed by a second round of SNP filtering (>5% missingness). Four samples at this stage failed quality control. Additionally, we removed SNPs that significantly deviated from the expected allele frequency assuming the Hardy-Weinberg equilibrium (P < 1 × 10⁻⁵). Last, SNPs were aligned with the TOPMed Hg38 Freeze 5 reference genome [17]. SNPs were checked for strand alignment, and variants with mismatching frequencies (R < 0.8) were removed in addition to AT-CG and INDEL (insertion-deletion) SNPs.

Sample Relatedness and Population Structure

By calculating identity by descent after pruning variants for linkage disequilibrium (R < 0.2) in PLINK, we estimated the genetic relatedness of our participants. One sample from each pair that shared alleles at a higher-than-expected frequency was removed (π > 0.185). In total, 4 sibling pairs were detected. To account for population structure during association analyses, the first 6 principal components were included, as determined by the elbow method [18] (Supplementary Figure 1). Outlying samples were removed if 2 SD from the mean principal component value, and 2 outliers were removed at this stage.

Imputation

Haplotypes were phased per the Eagle version 2 algorithm with the TOPMed Hg38 Freeze 5 reference genome, and Minimac4 imputed SNP data via the TOPMed server [19]. Variants were filtered by allele frequency concordance between the study data and the reference genome (R > 0.7) and refiltered for missingness and Hardy-Weinberg equilibrium. SNPs with a minor allele frequency <1% were excluded from the association analysis based on statistical power.

GWAS Modeling

To estimate if genotype was associated with RSV severity in our cohort, we fitted a generalized linear regression model of genotype vs RSV severity as determined by ReSVinet scores (0–20) using PLINK 2.0. Covariates were participant age, sex, preexisting conditions, including preterm birth, plus the first 6 principal components. The significance of the effect size of genotype on ReSVinet score was estimated by comparing the model slope with a null model using a likelihood ratio test under a chi-square distribution.

Functional Annotation of Candidate SNPs

SNPs identified in the GWAS of RSV severity that were significant above P < 1 × 10⁻⁴ were input into FUMA version 1.4.1 and SNPNexus variant annotation pipelines [20, 21]. Features included predicted variant deleteriousness and those previously reported in the GWAS catalog [22]. Chromatin sites and regulatory features of loci were examined with RegulomeDB and the Roadmap database [23, 24], and the Reactome database was used to examine pathway enrichment of SNPs, with the Fisher exact test to determine significance [21]. Variants with prior GWAS associations were queried in the GTEx version 8 database for tissue-specific eQTLs in whole blood and lung tissue [25].

Expression Quantitative Trait Loci Study

RNA Extraction and Transcriptomic Analysis

Blood samples were taken within 48 hours of admission, and RNA was extracted from whole blood with the QIAsymphony PAXgene Blood-RNA Kit (QIAGEN) according to the manufacturer's instructions and processed for Clariom GOScreen microarray (Thermo Fisher) after quality control. First-strand cDNA was synthesized by Poly-dT and random primers with a 5′-adaptor sequence and a 3′-adaptor for single-stranded cDNA, followed by low-cycle polymerase chain reaction amplification for labeling with the GeneChip Pico Reagent kit. In vitro transcription then amplified antisense mRNA (cRNA) for a second cDNA synthesis round, producing double-stranded cDNA. The targets were hybridized onto a single GO Screen plate after fragmentation, denaturation, and end labeling according to manufacturer instructions. A GeneChip Fluidics Station 450 and GeneChip Scanner 3000 7G stained and scanned single-sample cartridge arrays, while GeneTitan Multi-channel Instruments stained and imaged array plates.

eQTL Analysis

Genotypes from candidate SNPs (P < 1 × 10⁻⁴) were modeled with gene expression by a linear regression model via MatrixEQTL in R [26]. Covariates were age, sex, preexisting conditions (including premature birth), genetic principal components 1 to 6, plus transcriptome principal components 1 to 9, as determined by principal component analysis [18]. SNPs located within 1 Mb of the gene transcription start site were classified as cis-eQTLs. The false discovery rate was used to adjust for multiple testing; eQTLs of relevant variants were explored with R version 4.2; and statistical tests between groups were performed with the Wilcoxon test.

Code Availability

Full details of these scripts and methods can be found at https://github.com/mari-lynne/RSV_GWAS.

RESULTS

GWAS of RSV Severity in a Pediatric Cohort

Cohort Characteristics

Between 2017 and 2020, 261 infants infected with RSV were recruited across the United Kingdom, Netherlands, and Spain and enrolled for genomic and transcriptomic analyses following infection. The median age of our study cohort was ∼3 months (Table 1). Following diagnosis, participants were observed for clinical symptoms, such as temperature and wheeze, and need for mechanical ventilation, which were combined to give a composite ReSVinet score. Forty infants were recruited with preexisting conditions, including 14 preterm infants (<37 weeks).

When the genetic ancestry of our study cohort was evaluated, the majority of participants clustered within the European superpopulation (Supplementary Figure 1). To account for population structure, principal components 1 to 6 were included in the final regression model, and 2 participants with outlying distribution 2 SD from the mean principal components were removed. Four siblings were also removed to account for sample relatedness. Genotyping data were evaluated for SNP missingness and sample quality, for which an additional 4 samples were excluded. In total 251 participants were included from the initial 261 RSV cases used for GWAS analysis, with a median ReSVinet score of 8 (Table 1).

Host Genetic Factors Involved in RSV Infection Severity

To investigate genetic factors that may contribute to RSV severity, we conducted a GWAS using a linear model of ReSVinet, a clinical composite score that evaluates respiratory disease severity from a scale of 0 to 20. Age and previous comorbidities, including preterm birth, were included as covariates in the GWAS model, along with principal components 1 to 6. Following imputation, variants were filtered to include those with a minor allele frequency >0.01, which left 8 917 131 SNPs for downstream analyses. No variants were associated with RSV infection severity at the genome-wide significance level (P > 5 × 10⁻⁸). However, given the relatively small sample size, our study was underpowered to detect significant associations of moderate effect sizes (β = 0.2) per the standard genome-wide threshold [27].

When the Manhattan plot was subsequently examined, genetic loci with greater significance were discernible (Figure 1). Further examination of the distribution of P values under the observed and null hypothesis showed deviations from expected P values at those below P < 1 × 10⁻⁴ (Supplementary Figure 2). Given these data, we used a suggestive significance threshold of P < 1 × 10⁻⁴ for downstream analysis. In total, 816 variants associated with RSV severity were identified below the suggestive significance threshold of P < 1 × 10⁻⁴, and 75 variants had a P < 1 × 10⁻⁵.

Functional Annotation and Pathway Analysis of Candidate SNPs Associated With RSV Infection Severity

To identify variants that might be relevant in determining the severity of RSV infection, we explored genomic features and public GWAS data of candidate SNPs using SNPNexus and FUMA [20, 21]. Interestingly, the most significant SNPs surrounding the LSP1 gene locus on chromosome 11, region p15.5 (P < 1 × 10⁻⁶), were located upstream of the gene's promoter site (Supplementary Figure 4). Furthermore, LSP1 variants observed in our study were previously identified in GWASs of inflammatory bowel disease, ulcerative colitis, and white blood cell counts, suggesting a role for these SNPs in attenuating inflammation (Supplementary Table 1). Additionally, variant rs3843306 in chromosome 1, region p22.2, was associated with pulmonary function decline in those with asthma [28], while a chromosome 10 variant (10p12.33), rs2477664, in the MRC1 gene locus (macrophage mannose receptor) was reported in a GWAS of hematologic traits. Taken together, these results highlight potential candidate genes that may be relevant in the immune response during RSV infection.

To evaluate whether GWAS-associated genes were related to biologically relevant functions, we performed a pathway enrichment analysis of genes mapped to their closest variant within 10 kb of a transcription start or stop site [21]. An overall 113 GWAS-related genes were evaluated against genes annotated in the Reactome pathway database for overenrichment [29] (Supplementary Table 1). Immune system–related pathways, including antigen cross-presentation and butyrophilin family interactions (MRC1, PSMB1, BTNL8), were significantly enriched for GWAS-associated genes, while enriched disease pathways included genes involved in transmembrane support (ABCC6, PSMB1). The inactivation of CDC42 and RAC1 was also highlighted by the associations of SNPs in the SRGAP3 gene locus, a GTPase activator protein for RAC1.

eQTL Analysis of Candidate Variants Associated With RSV Infection Severity

We identified a significant proportion of our candidate variants within histone regions, as well as several SNPs that were within transcription factors (Supplementary Figure 3). To investigate the influence of these variants on gene transcription, we performed an eQTL analysis of candidate SNPs with gene expression measured in whole blood taken at clinical admission. In total, 816 SNPs were evaluated against approximately 1600 protein-coding genes measured on Affymetrix microarray chips. Cis-eQTLs included SNPs that were located within 1 Mb of a gene locus. Matrix eQTL was employed to test associations for gene transcription, with age, sex, and genetic principal components as model covariates.

After correction for multiple testing, we identified 2 variants in linkage disequilibrium (rs11071126, rs1501398) located upstream of the RAB27A gene that were significant eQTLs. These variants were associated with reduced RSV severity (GWAS, P < 1.39 × 10⁻⁵, β = −1.83) and a decrease in RAB27A expression in infected infants (adjusted P = .045, β = −0.29; Figure 2). As LSP1 gene expression was not significant for eQTL in the blood, we examined LSP1 variants (11p15.5) for associations with tissue-specific gene expression with the GTEx database [25]. Subsequently, the LSP1 variants were significantly associated with differences in LSP1 expression in lung tissue and whole blood (rs907612, P = 5.2 × 10⁻¹⁷, P = 2.36 × 10⁻⁹). SNP mapping revealed that these variants were located within regulatory regions of the gene locus, including transcription start sites (Supplementary Figure 4).

DISCUSSION

Here we present the first eQTL analysis of candidate SNPs associated with RSV severity, highlighting variants in the RAB27A gene locus that were associated with reduced RSV severity and lower gene expression. RAB27A is a small GTPase and a key regulator of neutrophil activation and exocytosis and is overexpressed during severe RSV [30, 31]. Moreover, in RAB7A/B knockout mouse models, neutrophil recruitment to the lung was significantly reduced [32]. Neutrophil infiltration in the lung is a hallmark of severe RSV infection; therefore, variants influencing these pathways are of particular interest. Additionally, genome annotation highlighted variants within the SRGAP3 locus associated with reduced RSV severity. SRGAP3 encodes another GTPase-activating protein linked to the inactivation of Rac1 and CDC42. Malhi et al recently observed that disrupting Rho GTPases and Rac1 activity through statins reduced RSV infection in vitro, highlighting potential therapeutic pathways to explore [33].

While no genome-wide significant hits were identified, we did observe several candidate variants associated with RSV severity that had been reported in GWASs of inflammatory traits, including variants in the MRC1, BTNL8, and LSP1 gene loci. Lymphocyte-specific protein 1 (LSP1) is primarily expressed in neutrophils, monocytes, and T cells, and its expression significantly increases during severe RSV pediatric infection, where it plays a crucial role in neutrophil activation and migration in the lung [31, 34]. During severe infections, increased recruitment of neutrophils in the respiratory tract can cause the release of proinflammatory cytokines, chemokines, and reactive oxygen species [35].

Subsequently, LSP1 gene knockouts in asthmatic mice reduced lung neutrophil infiltration, demonstrating how high levels of LSP1-induced neutrophil migration may be detrimental to disease progression [36]. In our cohort, we observed that GWAS variants of LSP1 were associated with differential gene expression in the lung and whole blood using data from the GTEx consortium [25]; however, it was not a significant eQTL in our data set. Larger sample sizes and potentially different tissue expression will therefore be needed to confirm these associations.

Pathway enrichment analysis of GWAS loci highlighted several immune-related genes, including MRC1, PSMB1, and BTNL8. The MRC1 gene encodes the mannose receptor C-type 1, which is involved in macrophage activation and class 1 antigen presentation [37]. In our study, the MRC1 variant rs2477664 was associated with reduced severity during RSV infection. Furthermore, examining single-cell sequencing data from lung tissue revealed that MRC1 and PSMB1 are highly expressed in alveolar macrophages, which are important in controlling RSV replication and progression through the respiratory tract [38, 39]. Therefore, it is possible that early antigen presentation might be an important factor in controlling RSV disease. Interestingly, Pasanen et al identified genes in overlapping pathways, including antigen presentation, and highlighted genes also present in our study, such as DSCAM and TRABD2B, which have yet to be explored in the context of RSV infection [11].

The main limitation of this study was the small sample size relative for a GWAS and the population's homogeneity. Given the high number of comparisons made across hundreds of thousands of SNPs, large sample sizes are typically needed to detect robust associations with disease. Consequently, examination of these findings in new cohorts will be critical to confirming these results. Two other GWAS's of RSV severity had been conducted at the time of writing, with bronchiolitis as the diagnostic criterion for severe RSV in one study and prolonged RSV infection in the other [11, 40]. While neither of these studies detected any genome-wide significant results, they may be useful for further meta-analyses. Future research may continue to investigate the role of GWAS-associated genes in the host-immune response to RSV, particularly those linked to neutrophil pathways and antigen presentation, shedding light on potential therapeutic avenues.

Supplementary Data

Supplementary materials are available at The Journal of Infectious Diseases online (http://jid.oxfordjournals.org/). Supplementary materials consist of data provided by the author that are published to benefit the reader. The posted materials are not copyedited. The contents of all supplementary data are the sole responsibility of the authors. Questions or messages regarding errors should be addressed to the author.

Notes

PROMISE investigators. Harish Nair (University of Edinburgh), Hanna Nohynek (THL, Finland), Anne Teirlinck (RIVM, The Netherlands), Louis Bont (University Medical Center Utrecht), Peter Openshaw (Imperial College, London), Andrew Pollard (University of Oxford), Veena Kumar (Novavax), Tin Tin Htar (Pfizer), Charlotte Vernhes and Rolf Kramer (Sanofi Pasteur), Gael Dos Santos (GlaxoSmithKline), Jeroen Aerssens (Janssen).

Acknowledgments. This work was supported by the Respiratory Syncytial Virus Consortium in Europe (RESCEU) project and Preparing for RSV Immunisation and Surveillance in Europe (PROMISE).We would like to acknowledge the following PROMISE investigators for their contributions to this study.; Hanneke van Zoggel (University Medical Center Utrecht), Sarah Hak (University Medical Center Utrecht), Rob van Binnendijk (National Institute for Public Health and the Environment), Berta Gumí Audenis (Teamit Research S.L.), Maica Llavero (Teamit Research S.L.), Cornelia Gottschick (Martin Luther Universitaet Halle-Wittenberg), Rafael Mikolajczyk (Martin Luther Universitaet Halle-Wittenberg), Annick Moureau (Sanofi Pasteur SA), Donghui Zhang (Sanofi Pasteur SA), Jim Janimak (GlaxoSmithKline Biologicals SA), Martin Ota (GlaxoSmithKline Biologicals SA), Gabriela Ispas (Janssen Pharmaceutica NV), Paul Peeters (Janssen Pharmaceutica NV), Rosemary Bastian Arangassery (Janssen Pharmaceutica NV), Hadi Beyhaghi (Novavax INC), Beate Schmoele-Thoma (Pfizer Limited), Kena Swanson (Pfizer Limited), Kimberly Shea (Pfizer Limited), and Reiko Sato (Pfizer Limited). Also, we thank Larry Anderson (Emory University) for reviewing the draft manuscript.

Disclaimer. This article has not been submitted elsewhere, including preprint servers, nor have the data been presented publicly.

Financial support. This work was supported by the RESCEU project and PROMISE (Preparing for RSV Immunisation and Surveillance in Europe); the Innovative Medicines Initiative 2 Joint Undertaking (grants 116019 and 101034339 to RESCEU and PROMISE, respectively); and the European Union’s Horizon 2020 Research and Innovation Programme and the European Federation of Pharmaceutical Industries and Associations (to the Innovative Medicines Initiative 2 Joint Undertaking).

Supplement sponsorship. This article appears as part of the supplement “Preparing Europe for Introduction of Immunization Against RSV: Bridging the Evidence and Policy Gap.”

References

Author notes