Genome-wide association study of long COVID

Abstract

Infections can lead to persistent symptoms and diseases such as shingles after varicella zoster or rheumatic fever after streptococcal infections. Similarly, severe acute respiratory syndrome coronavirus 2 (SARS‑CoV‑2) infection can result in long coronavirus disease (COVID), typically manifesting as fatigue, pulmonary symptoms and cognitive dysfunction. The biological mechanisms behind long COVID remain unclear. We performed a genome-wide association study for long COVID including up to 6,450 long COVID cases and 1,093,995 population controls from 24 studies across 16 countries. We discovered an association of FOXP4 with long COVID, independent of its previously identified association with severe COVID-19. The signal was replicated in 9,500 long COVID cases and 798,835 population controls. Given the transcription factor FOXP4's role in lung physiology and pathology, our findings highlight the importance of lung function in the pathophysiology of long COVID.

Fig. 1. Geographic overview of studies contributing to the Long COVID HGI.

The 24 studies contributing to the Long COVID HGI data freeze 4 served as the discovery cohorts for the GWAS meta-analyses. Each color represents a meta-analysis with specific case and control definitions. Strict case definition, long COVID after test-verified SARS-CoV-2 infection; broad case definition, long COVID after any SARS-CoV-2 infection; strict control definition, individuals that had SARS-CoV-2 but did not develop long COVID; broad control definition, population control, that is, all individuals in each study that did not meet the long COVID criteria. Effective sample sizes are shown as the size of each diamond shape, and locations of sample collection in (from left to right) North America, Europe, Middle East and Asia. For more detailed sample sizes, see Supplementary Table 11.

Fig. 2. Meta-analysis of 11 GWAS studies of long COVID shows an association at the FOXP4 locus.

a, Manhattan plot of long COVID after test-verified SARS-CoV-2 infection (strict case definition, n = 3,018) compared to all other individuals in each dataset (population controls, broad control definition, n = 994,582). A genome-wide significant association with long COVID was found in the chromosome 6, upstream of the FOXP4 gene (chr6: 41,515,652 G:C, GRCh38, rs9367106, as the lead variant; P = 1.76 × 10⁻¹⁰, Bonferroni P = 7.06 × 10⁻¹⁰, increased risk with the C allele, OR = 1.63, 95% CI = 1.40–1.89). Horizontal lines indicate genome-wide significance thresholds for IVW meta-analysis before (P < 5 × 10⁻⁸, dashed line) and after (1.25 × 10⁻⁸) Bonferroni correction over the four long COVID meta-analyses (INCMNSZ = MexGen-COVID Initiative). b, Chromosome 6 lead variant across the contributing studies and ancestries in GWAS meta-analyses of long COVID with strict case definition and broad control definition. Lead variant rs9367106 (solid line) and if missing, imputed by the variant with the highest LD with the lead variant for illustrative purpose, that is, rs12660421 (r = 0.98 in European in 1,000 G + HGDP samples, dotted lines). For the imputed variants, β was weighted by multiplying by the LD correlation coefficient (r = 0.98). Centre, OR; error bar, 95% CI. Genetic ancestries marked by colors. MAF varies across ancestries, ranging from 1% to 34% (Supplementary Fig. 4). AFR, African; AMR, Admixed American; EAS, East Asian; EUR, European; UKBB, UK Biobank. (Results for the other three GWAS meta-analyses in Supplementary Figs. 2 and 3a–c).

Fig. 3. The chromosome 6 region (chr6: 41,490,001–41,560,000 (70 kb); FOXP4 locus) in the long COVID GWAS meta-analysis.

Long COVID meta-analysis with strict case (n = 3,018) and broad control (n = 994,582) definition (Fig. 2). X axis shows the position on chromosome 6 (GRCh38). The long COVID lead variant (rs9367106) is depicted with a triangle in each plot. a, Locus zoom plot with each variant colored by effective sample size and showing statistical significance (IVW GWAS meta-analysis −log₁₀ P value) on y axis. b, Each variant colored by statistical significance and showing effect sizes (center, coefficients; error bar, 95% CI on y axis). c, Each variant colored by ancestry and showing LD correlation coefficient (r) with the long COVID lead variant on y axis. d, Ensembl genes in the region (FOXP4 not fully shown; www.ensembl.org).

Fig. 4. FOXP4 expression in the lung.

a, The lead variant rs9367106 was not found in the GTEx dataset, but a proxy variant (rs12660421, chr6: 41,520,640) in high LD (r² = 0.97, rs12660421-A allele is correlated with the long COVID risk allele rs9367106-C) showed a significant eQTL after multiple testing correction, increasing FOXP4 expression in the lung (P = 5.3 × 10⁻⁹, NES; expression with GA genotype compared to expression with GG, normalized to 0) = 0.56; GTEx V8 lung samples with GG genotype, n = 483, GA genotype, n = 32; https://gtexportal.org/home/snp/rs12660421). For other tissues, see multitissue eQTL plot in Supplementary Fig. 6. b, Colocalization analysis using eQTL data from GTEx v8 tissue type and long COVID GWAS meta-analysis association data (Supplementary Note). Plots illustrate −log₁₀ P value for long COVID (x axis) and for FOXP4 expression in the lung (y axis), regional association of the FOXP4 locus variants with long COVID (top right) and regional association of the FOXP4 variants with RNA expression measured in the lung in GTEx (bottom right). Variants are colored by 1000 Genomes European-ancestry LD r² with the lead variant (rs12660421) for FOXP4 expression in lung tissue (the most significant long COVID variant overlapping the GTEx v8 dataset (rs9381074) also annotated). c, Human Protein Atlas RNA single-cell type tissue cluster data (transcript expression levels summarized per gene and cluster) of lung (GSE130148) showing FOXP4 expression in unaffected individuals. The values were visualized using log₁₀ (pTPM + 1) values. Each annotation is taken from the clustering results performed in the Human Protein Atlas. pTPM, protein transcripts per million.

Fig. 5. Genetic correlations and MR causal estimates between long COVID and potential risk factors, biomarkers and diseases.

a,b, LD score regression (a, LDSC, top; Supplementary Table 27) and IVW MR (b, fixed-effects model,bottom; Supplementary Table 29 and Supplementary Data) were used for calculating two-sided P values. The size of each colored square corresponds to statistical significance (***P < 0.0001, full-sized square; **P < 0.01, full-sized square; *P < 0.05, full-sized square; P < 0.1, large square; P < 0.5, medium square and P > 0.5, small square; not corrected for multiple comparisons). A full list of traits is provided in Supplementary Table 26. For sample sizes in each long COVID GWAS meta-analysis using strict (S) or broad (B) case and control definitions, see Supplementary Table 11. c, MR scatter plot with effect sizes (β ± s.e.) of each variant on COVID-19 susceptibility (reported SARS-CoV-2 infection) as exposure and long COVID (strict case, broad control definition) as outcome (P (IVW, fixed effects) = 1.8 × 10⁻⁷, pleiotropy P = 0.47; Supplementary Table 30). d, Similarly, MR with COVID-19 hospitalization as exposure and long COVID as outcome (P (IVW fixed effects) = 4.8 × 10⁻⁸, pleiotropy P = 0.83; Supplementary Table 30). e, Analysis of shared and unique effects between SARS-CoV-2 infection susceptibility and long COVID using a Bayesian mixture model showed ABO and 3p21.31 rs73062389 as having shared effects (posterior probability > 0.99). FOXP4 variant association was discovered in the long COVID meta-analyses but showed also an effect on the susceptibility of the initial infection, though smaller than on long COVID (Supplementary Table 34). (Effects shown as β, error bars represent 95% confidence intervals.) f, Similarly, analysis of shared and unique effects between COVID-19 severity and long COVID using a Bayesian mixture model showed FOXP4 variant with a joint effect (posterior probability > 0.9), differing from the other severity variants due to its larger effect on long COVID (Supplementary Table 35). BMI, body mass index; CRP, C-reactive protein; eGFR, estimated glomerular filtration rate; ADHD, attention-deficit hyperactivity disorder.