No Genetic Causal Relationship between COVID-19 and Herpes Zoster: A Bidirectional Mendelian Randomization Analysis

Introduction

The COVID-19 pandemic, caused by the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2), has led to unprecedented global health challenges. Beyond its immediate respiratory implications, COVID-19 has been associated with a wide array of systemic complications, including cardiovascular, neurological, autoimmune, and dermatological conditions.^1–6 We recently reported an atypical herpes zoster (HZ) case presenting with urinary tract irritation as the initial symptom in a COVID-19 patient,⁷ which raises our concern of whether any potential interplay exists between HZ and COVID-19.

HZ, commonly known as shingles, is caused by the reactivation of the varicella-zoster virus (VZV) in individuals with a history of chickenpox. HZ typically manifests as a painful, vesicular rash along a dermatome and can lead to severe complications such as postherpetic neuralgia, especially in older adults and immunocompromised individuals.^8,9 The reactivation of VZV is thought to be linked to a decline in cell-mediated immunity, which can be influenced by various stressors, including infections.¹⁰ Given the immunomodulatory effects of SARS-CoV-2 and the immune dysregulation observed in COVID-19 patients,¹ there is a plausible biological basis for investigating a potential causal relationship between COVID-19 and herpes zoster reactivation. Recent studies have highlighted the increased incidence of HZ following COVID-19.^11–15 Additionally, a prospective cohort study found that younger individuals (median age: 42.5 years) are more likely to develop HZ after COVID-19 infection, while females show a higher propensity for HZ following COVID-19 vaccination.¹⁶ These observations suggest a possible association between COVID-19 and the onset of HZ.

Mendelian randomization (MR) offers a powerful epidemiological approach to assess causal relationships between exposures and outcomes using genetic variants as instrumental variables.¹⁷ This method mitigates confounding and reverse causation biases, providing more robust causal inferences compared to traditional observational studies.^18,19 By leveraging genetic variants associated with COVID-19 and HZ, MR can help elucidate the potential genetic causal relationship between these two diseases.

In this study, we aim to analyze the causal effect between COVID-19 and HZ using a bidirectional MR analysis. We utilized publicly available genome-wide association study (GWAS) summary statistics for COVID-19 and HZ to identify relevant genetic instruments. Through rigorous MR analysis, we seek to clarify whether COVID-19 serves as a causal trigger for HZ, or vice versa, in hope to provide further insights to the understanding of post-infectious sequelae.

Materials and Methods

Study Design

We conducted a bidirectional Mendelian randomization (MR) analysis to investigate the causal relationship between COVID-19 and HZ as shown in Figure 1. MR utilizes genetic variants associated with an exposure as instrumental variables (IVs) to infer causality regarding an outcome. The IVs were extracted strictly following the three assumptions: (1) IVs are related to the exposure; (2) IVs are independent of confounders; and (3) IVs affect the outcome only through the exposures.

Figure 1 Study design of this bidirectional MR analysis.

Data Sources

The MR analysis aimed to elucidate the potential causal relationship between COVID-19 and HZ, utilizing genetic instruments extracted from two major databases as listed in Table 1. The study populations were predominantly of European ancestry to ensure genetic homogeneity and minimize population stratification. The inclusion and exclusion criteria for our study aligned with those of the respective datasets, while only the data from European populations of the COVID-19 host Genetics Initiative extracted to maintain consistency in genetic background.

Table 1 The Summary Data of All GWAS Used in This Study

Covid-19

The genetic variants associated with COVID-19 were sourced from the largest GWAS for COVID-19, specifically the Round 7 release of the COVID19-hg GWAS meta-analyses by the COVID-19 host Genetics Initiative (https://www.covid19hg.org/results/r7/).^20–22 This dataset, a meta-analysis of up to 219,692 cases and over 3 million controls, encompasses 82 studies from 35 countries, including 36 studies of individuals with non-European ancestry. Three phenotypes were analyzed: (1) COVID-19 infection, representing general susceptibility (Susceptibility, 219,692 cases); (2) Hospitalized COVID-19 cases, indicating moderate to severe disease manifestations (Hospitalization, 49,033 cases); and (3) Very severe respiratory confirmed COVID-19 cases, highlighting the most critical disease outcomes (Severity, 21,193 cases).

Hz

As for HZ, the trait of “herpes zoster” was searched in the IEU GWAS database (https://gwas.mrcieu.ac.uk), and the associated genetic variants in the dataset “finn-b-AB1_ZOSTER” were obtained. This dataset, as part of FinnGen project with 500,000 Finnish biobank donors to understand the genetic basis of diseases (https://www.finngen.fi/en), includes 2080 hZ cases and 211,856 controls, with totaling 16,380,433 SNPs, all of European ancestry.

Instrument Selection and Linkage Disequilibrium (LD) Clumping

The selection of Single Nucleotide Polymorphisms (SNPs) as instruments for the exposure followed stringent criteria to ensure the relevance and robustness of the MR analysis. SNPs were filtered based on a significance threshold of p < 5×10⁻⁸ for the COVID-19 categories; and for HZ, a relaxed threshold of p < 1×10⁻⁵ was applied to ensure a sufficient number of SNPs for robust genetic analyses. Besides, LD clumping was performed to address the issue of linkage disequilibrium that could bias the MR estimates. The LD clumping parameters were set with a window of 10,000 base pairs, an r² cutoff of 0.001 to ensure the independence of SNPs, and significance levels for index SNPs and secondary SNPs both set at 1, utilizing data from the European population to maintain consistency and relevance to the genetic background of the study cohort.

MR Analysis

The “TwoSampleMR” package (https://github.com/MRCIEU/TwoSampleMR) in R (RStudio 2024.04.2 Build 764 with R 4.4.1) was utilized for conducting the MR analysis, employing Inverse-variance weighted (IVW) method as the primary, and MR-Egger, Weighted median, Weighted mode and Simple mode methods as the supplementary. If no horizontal pleiotropy existed, fixed effects IVW was used; otherwise, random effects IVW was employed. The relationship between COVID-19 genetic risk and HZ incidence was illustrated through scatter plots, while the individual and combined effects of SNPs were depicted in forest plots.

Results were expressed as beta coefficients with standard errors (β ± SE), and a P-value less than 0.05 was considered indicative of a statistically significant difference.

Results

Instrument Extraction

To access the potential causal effect of COVID-19 on HZ, we extracted three categories of COVID-19 SNPs, susceptibility, hospitalization, and severity, as IVs from the round 7 release of the COVID19-hg GWAS meta-analyses.²⁰ After clumping, we retrieved 29 SNPs for COVID-19 susceptibility, 33 SNPs for COVID-19 hospitalization, and 15 SNPs for COVID-19 severity (Table S1-3).

Conversely, to evaluate the causal effect of HZ on COVID-19 outcomes, we obtained HZ-associated SNPs from the IEU Open GWAS database using the phenotype identifier “finn-b-AB1_ZOSTER”. Following clumping procedures, 19 SNPs were retained as IVs for HZ (Table S4).

No Significant Causal Effect between COVID-19 and HZ was Observed

We conducted a bidirectional MR analysis to investigate the potential causal relationship between COVID-19 and HZ, assessing the effects of COVID-19 susceptibility, hospitalization, and severity on HZ and vice versa. We utilized the inverse variance weighted (IVW) method as the primary analysis, supplemented by MR-Egger, weighted median, simple mode, and weighted mode methods, to strengthen our findings (Figure 2 and Table S5).

Figure 2 Bidirectional MR analysis results indicated no causal effect between COVID-19 and HZ. (A-C) The forward MR analysis results of the causal effect of COVID-19 ((A) susceptibility; (B) hospitalization, (C) severity) on HZ. (D-F) The reverse MR analysis results of the causal effect of HZ on COVID-19 susceptibility (D), hospitalization (E) and severity (F). The results are shown as the effect size and the standard error (β±SE).

In the forward MR analysis, where COVID-19 traits were considered as exposure and HZ as the outcome, the IVW method showed no significant causal effect: COVID-19 susceptibility on HZ (β = −0.053, se = 0.182, p = 0.77), COVID-19 hospitalization on HZ (β = 0.060, se = 0.069, p = 0.38), and COVID-19 severity on HZ (β = 0.015, se = 0.048, p = 0.75), which were consistent across other MR methods (Figure 2A–C and Table S5). Scatter plots visualizing SNP effects for COVID-19 traits against those for HZ (Figure 3A–C) showed slopes close to zero, indicating minimal or no causal association. Forest plots of individual SNPs (Figure 3D–F) reinforced the absence of a significant association, as indicated by combined effects from IVW and MR-Egger.

Figure 3 Forward MR analysis of the causal effect of COVID-19 on HZ. (A-C) The scatter plots illustrating the relationship between SNP effects on COVID-19 ((A) susceptibility; (B) hospitalization, (C) severity) and SNP effects on HZ. (D-F) The forest plots displaying the Wald ratio for single SNPs of COVID-19 susceptibility (D), hospitalization (E) and severity (F), as well as their combined effect calculated by MR-Egger and IVW methods.

In the reverse MR analysis, where HZ was treated as the exposure and COVID-19 traits as the outcomes, the IVW results similarly showed no significant causal effect: HZ on COVID-19 susceptibility (β = 0.006, se = 0.006, p = 0.33), HZ on COVID-19 hospitalization (β = −0.012, se = 0.012, p = 0.32), and HZ on COVID-19 severity (β = −0.015, se = 0.020, p = 0.46) (Figure 2D–F and Table S5). Notably, the MR-Egger method detected a significant causal effect of HZ on COVID-19 hospitalization (β = −0.042, se = 0.017, p = 0.02). However, given the low β value and near-zero slopes in scatter plots (Figure 4A–C), along with forest plot results (Figure 4D–F), this effect is likely negligible, and overall evidence does not support a causal relationship between HZ and COVID-19 outcomes.

Figure 4 Reverse MR analysis of the causal effect of HZ on COVID-19. (A-C) The scatter plots illustrating the relationship between SNP effects on HZ and SNP effects on COVID-19 ((A) susceptibility; (B) hospitalization, (C) severity). (D-F) The forest plots displaying the Wald ratio for single SNPs of HZ on COVID-19 susceptibility (D), hospitalization (E) and severity (F), as well as their combined effect calculated by MR-Egger and IVW methods.

Sensitivity Analysis

We evaluated the robustness of the MR findings through sensitivity analyses. The heterogeneity among the IVs were assessed using Cochran’s Q tests in both IVW and MR-Egger models, no significant heterogeneity for either forward or reverse analyses (Table S6). Visual inspection of funnel plots indicated symmetry in most analyses (Figure 5A–C and F), suggesting minimal directional pleiotropy. However, asymmetry was observed in plots for HZ on COVID-19 susceptibility and hospitalization (Figure 5D and E), indicating possible pleiotropic effects that might influence the results for HZ on COVID-19 susceptibility and hospitalization.

Figure 5 Funnel plots for MR analyses examining the bidirectional association between COVID-19 and HZ. (A-C) The forward MR analysis, where COVID-19 ((A) susceptibility; (B) hospitalization, C. severity) are considered as exposures and HZ as the outcome. (D-F) The reverse MR analysis, where HZ is the exposure and COVID-19 susceptibility (D), hospitalization (E) and severity (F) are the outcomes.

Then, the MR-radial plots were used to identify the potential outliers among IVs with MR-Egger and IVW methods. As shown in Figure 6, no outliers were observed in COVID-19 IVs on HZ (Figure 6A–C) or in HZ IVs on COVID-19 susceptibility (Figure 6D). One SNP (rs2523580) was identified as an outlier in HZ IVs on COVID-19 hospitalization and severity (Figure 6E and F).

Figure 6 MR-radial plots for the bidirectional MR analyses. (A-C) The forward MR analysis, where COVID-19 ((A) susceptibility; (B) hospitalization, (C) severity) are considered as exposures and HZ as the outcome. (D-F) The reverse MR analysis, where HZ is the exposure and COVID-19 susceptibility (D), hospitalization (E) and severity (F) are the outcomes. The blue line represents the inverse variance weighted (IVW) method, while the Orange line represents the MR-Egger method. Outlier SNPs are highlighted with circles, as shown in panels E and F.

Horizontal pleiotropy was further assessed using the MR-Egger intercepts and MR-PRESSO global tests, and the results are shown in Table S6. In the forward MR analysis, MR-Egger regression intercepts for COVID-19 susceptibility (intercept = 0.004, p = 0.83), hospitalization (intercept = 0.005, p = 0.71) and severity (intercept = 0.001, p = 0.96) on HZ was not significant, consistent with the MR-PRESSO global test results (COVID-19 susceptibility, RSS_obs = 11.093, p = 0.80; hospitalization, RSS_obs = 37.588, p = 0.32; severity, RSS_obs = 31.917, p = 0.39).

Of note, in the reverse MR analysis, minimal pleiotropy was identified by significant MR-Egger intercepts for COVID-19 susceptibility (intercept = 0.005, p = 0.04) and hospitalization (intercept = 0.012, p = 0.03), but not for severity (intercept = 0.001, p = 0.90). However, the MR-PRESSO global tests indicate no significant horizontal pleiotropy (COVID-19 susceptibility, RSS_obs = 23.683, p = 0.28; hospitalization, RSS_obs = 21.865, p = 0.37; severity, RSS_obs = 25.374, p = 0.19), suggesting the pleiotropy effects were limited.

Potential outliers among the SNPs were further examined through the leave-one-out analysis, re-estimating causal effects by sequentially excluding each SNP. No significant changes in estimates were observed when SNPs were removed individually (Figure 7), suggesting robustness of the MR results to outliers.

Figure 7 Leave-one-out sensitivity analysis. (A-C) Forward MR analysis with COVID-19 ((A) susceptibility; (B) hospitalization, C. severity) as the exposure and HZ as the outcome. (D-F) Reverse MR analysis with HZ as the exposure and COVID-19 (D) susceptibility, (E) hospitalization and (F) severity) as the outcome. The results are the causal estimates and 95% confidence intervals (CIs) of the MR analysis after omitting each individual SNP.

Discussion

This MR study aimed to elucidate the potential causal relationship between COVID-19 and HZ. Despite the well-documented immunological impact of COVID-19 and the clinical observations suggesting a higher incidence of HZ among COVID-19 patients, our analysis found no significant genetic causal effect. This lack of association was consistent across multiple MR methods, including MR Egger, weighted median, inverse variance weighted (IVW), simple mode, and weighted mode approaches.

The absence of a causal link in our MR analysis implies that the genetic predisposition to COVID-19 does not directly influence the risk of HZ reactivation or vice versa, suggesting that the observed clinical correlation between COVID-19 and herpes zoster may be due to confounding factors rather than a direct genetic causal pathway. This is a critical insight, as it suggests that other factors, such as the overall immune status of patients, concurrent treatments, or stress levels associated with the pandemic, might play more substantial roles in the increased incidence of HZ observed clinically.^11–15 For instance, the immunosuppressive environment induced by severe COVID-19 or its treatments might provide a conducive context for varicella-zoster virus (VZV) reactivation without there being a direct genetic causative pathway. Severe COVID-19 is associated with lymphopenia and T-cell exhaustion, reducing CD8+ T-cell responses critical for maintaining VZV latency.²³ Additionally, treatments such as corticosteroids and immunosuppressive therapies, used to manage severe cases, further impair T-cell function,²⁴ increasing HZ risk. Pandemic-related stress may also exacerbate immune suppression, facilitating VZV reactivation.²⁵

Previous studies, including a retrospective cohort study by Chen et al indicated a higher risk of herpes zoster following COVID-19, but it could not eliminate potential confounders.¹⁵ Additionally, Diez-Domingo et al reviewed 27 reported cases of HZ following COVID-19, and emphasized the need for awareness among practitioners about the potential increased risk of HZ during the COVID-19 pandemic and recommended considering preventive measures.²⁶ Algaadi et al reviewed available literature on HZ in COVID-19 patients, highlighting the need for larger and more comprehensive studies to confirm a causal relationship between COVID-19 and HZ.²⁷ Our findings provide important context for these studies, suggesting that the link between COVID-19 and HZ may not be genetically mediated but instead may be influenced by external or environmental factors.

The strengths of our study include the use of large-scale GWAS datasets and the application of various MR methods to ensure the robustness of our findings. However, several limitations must be considered. First, our analysis was restricted to populations of predominantly European ancestry, which might limit the generalizability of our findings to other ethnic groups. Second, the genetic instruments used, while robust, may not capture all aspects of COVID-19 traits, potentially missing more nuanced genetic influences. Lastly, the possibility of residual confounding cannot be entirely ruled out, although MR is designed to minimize such biases.

The clinical implication of our findings is that preventive and therapeutic strategies for HZ in the context of COVID-19 should consider factors beyond genetic predisposition. Healthcare providers should remain vigilant for HZ in COVID-19 patients, especially those with known immunosuppressive conditions or treatments. Future research should focus on exploring the non-genetic factors contributing to the increased herpes zoster incidence observed during the pandemic. Larger, more diverse population studies and investigations into the immunological mechanisms linking COVID-19 and herpes zoster are warranted to fully understand the interplay between these conditions.

Conclusion

Our MR analysis found no evidence of a causal relationship between COVID-19 and HZ, suggesting that the increased incidence of HZ observed among COVID-19 patients in clinical settings may result from other factors rather than a direct genetic link. These findings underscore the importance of a comprehensive approach to understanding and managing the long-term effects of COVID-19 on health, considering both genetic and non-genetic influences.