Genetic Evidence Linking Immune Cell Subsets to Psoriatic Arthritis Susceptibility: A Mendelian Randomization Study

Introduction

Psoriatic arthritis (PsA) is an immune-mediated inflammatory arthropathy that manifests as seronegative spondyloarthritis in patients with psoriasis.¹ Its clinical presentation encompasses psoriatic rash, joint swelling, pain, stiffness, and impaired joint mobility, accompanied, in certain cases, by spondylitis and sacroiliac arthritis,² which contributes to a substantial burden on both families and healthcare systems. One of the factors making it challenging to identify and treat psoriatic arthritis in its early stages is its unclear etiology.³ Consequently, there is a pressing need to identify novel risk factors for PsA to facilitate the development of preventative and therapeutic strategies, thereby alleviating the future burden on the healthcare system.

Within the realm of immunology, leukocytes, which comprise lymphocytes and diverse phagocytic cells, are instrumental in orchestrating immune responses. Furthermore, studies have demonstrated that immune cells significantly influence multiple processes in bone physiology, including bone regeneration, osteoclastogenesis, osteoblast activity, bone density modulation, and other pertinent functions associated with skeletal homeostasis.^4,5 Alterations in immune cells have been linked to the risk, severity, and progression of PsA. For instance, one study identified monocytes, dendritic cells, memory T regulatory cells, and T-helper 17 cells as being associated with increased disease activity and longer PsA duration.⁶ Evidence also suggests that disrupted homeostasis of innate lymphoid cells may confer a potential risk of PsA.⁷ Therefore, using MR analysis, we examined the complicated causative link as well as the reverse causal association between 731 immune cell types and immunocyte phenotypes in PsA.⁸ The insights generated by this study may inform future research into PsA mechanisms and the development of clinical diagnostics and therapies.

Recent progress in computational analytic methods has introduced new approaches for examining relationships between distinct disease phenotypes.^9,10 These techniques offer valuable opportunities to clarify the mechanistic basis of PsA and to uncover candidate targets for therapeutic intervention. Although randomised controlled trials are commonly recognised as the most effective method for determining causal associations in epidemiological studies, their practical execution is often constrained by ethical considerations and substantial financial outlays.¹¹ Studies based on observation face inherent challenges from unmeasured confounders and bidirectional causality, which undermine confident causal conclusions.¹² Mendelian randomization (MR) addresses these challenges through the use of genetic variants as instrumental proxies, enabling more reliable causal estimation.^13–15 Jia et al (2025) previously conducted a bidirectional MR study on 731 immune phenotypes in PsA,⁸ and they concluded that there was no reverse causal relationship⁸—this finding may have overlooked the immunorenewal phenomenon that has been confirmed in clinical PsA. Our study expands and deepens this previous work in two fundamental aspects: Firstly, by providing a more detailed description of the reverse causal path, our research results are coordinated with the clinical evidence of systemic immune changes (such as the alterations of Th17 and regulatory T cells).^16–18 Secondly, by isolating the unified biological characteristics of dysregulation of B cell co-stimulatory markers—especially CD25 on multiple memory B cell subpopulations. In contrast to prior work, we pay particular attention to the biological plausibility of B-cell co-stimulatory signals and offer a more circumspect interpretation of the reverse MR findings in the context of disease-driven immune remodeling. Our objective was to pinpoint immune markers that might serve as targets for early detection or therapeutic development.

Method

Study Design

We employed bidirectional two-sample MR analysis to investigate the causal relationship between 731 immunocyte phenotypes and PsA. First, the direct effects of each immunocyte phenotype on the occurrence of PsA were estimated using MR analyses (Figure 1). Three main hypotheses were tested: (1) The IVs displayed a strong correlation with immunocyte phenotypes, (2) the IVs remained unaffected by confounding factors related to the association between immunocyte phenotypes and PsA, and (3) the IVs had no effect on immunocyte phenotypes through pathways other than immunocyte phenotypes. Subsequently, we performed reverse MR analyses treating PsA as the exposure and immune cell phenotypes as outcomes, thereby testing whether disease liability influences immune composition. All analyses adhered to the STROBE-MR reporting guidelines. Effect estimates are expressed as odds ratios (ORs) for PsA per unit change in immunocyte phenotypes.

Figure 1 Schematic illustration of the core assumptions underpinning Mendelian randomisation. Solid arrows (→) denote the hypothesised causal pathway from instrumental variables (IVs; SNP1–SNPn) through the exposure (immunocyte phenotype) to the outcome (psoriatic arthritis). Dashed arrows (--→) represent potential violations of the Independence and exclusion restriction assumptions via confounders (nutritional status, body mass index, obesity, smoking, and excessive alcohol consumption). Red cross symbols (×) indicate pathways that must be absent for valid causal inference: specifically, the IVs should not be associated with confounders (independence assumption) and should not affect the outcome through any pathway other than the exposure (exclusion restriction assumption).

Data Source

Exposure data for 731 immune cell phenotypes were obtained from a large-scale genome-wide association study (GWAS; GWAS Catalog accession IDs: GCST90001391–GCST90002121) conducted in 3,757 individuals of European ancestry, enabling the construction of a standardized, high-resolution atlas of the human immune system. For the outcome analysis, we selected the PsA dataset from the IEU OpenGWAS database (ID: IEU-B-5116, https://gwas.mrcieu.ac.uk/datasets), comprising 5,065 cases and 21,286 controls—all of European ancestry—yielding a total sample size of N = 26,351.¹⁹ This dataset was prioritized following systematic evaluation of key methodological criteria, including sample size, single-nucleotide polymorphism (SNP) density, and ancestral homogeneity, to maximize statistical power and mitigate potential confounding due to population stratification.

Genetic Data Processing

We employed the R package TwoSampleMR (version 0.5.6 in R 4.3.0) to identify valid genetic instruments from exposure summary statistics. To ensure independence among selected variants and minimize confounding from linkage disequilibrium, we applied a clumping procedure using an r² threshold of 0.001 within a 10 Mb window, referencing the 1000 Genomes European ancestry panel.²⁰ Given the limited sample size of the immune cell GWAS, we adopted a relaxed genome-wide significance threshold (P < 5×10^–5) for instrument selection, consistent with recent MR studies examining immune traits. However, we acknowledge that this relaxed threshold introduces a methodological trade-off: while it increases the number of available instruments and improves statistical power, it may also admit variants with weaker exposure associations, elevating the potential for false-positive findings. To mitigate this risk, we applied stringent post-selection quality control, including F-statistic filtering and Steiger directionality testing, and interpret all findings in the context of this limitation. To exclude weak instruments susceptible to bias, we calculated F-statistics using the formula F = β²/SE², retaining only variants with F > 10.²¹ This conventional threshold was selected because F-statistics below 10 are empirically associated with meaningful weak instrument bias in MR settings. Across all retained instruments, F-statistics ranged from 29.85 to 5062, indicating that the majority of selected variants satisfied instrument strength requirements. Following instrument selection, we harmonized exposure and outcome datasets by aligning effect alleles to a common reference strand.²² Ambiguous palindromic SNPs lacking clear strand orientation were removed, while those with unambiguous allele configurations were appropriately oriented to ensure directional consistency.

MR Analyses Between 731 Immunocyte Phenotypes and PsA

To investigate putative causal relationships, we conducted MR analyses using multiple complementary approaches: IVW analysis under both fixed-effects and random-effects frameworks, weighted median estimation, MR-Egger regression, and weighted mode analysis. The IVW method served as our primary analytical framework, with Cochran’s Q statistic used to evaluate heterogeneity across genetic instruments. Given the large number of independent tests (n = 731), we applied Benjamini-Hochberg correction to control the false discovery rate (FDR). Associations surviving FDR correction (FDR < 0.05) were considered robust findings. Associations achieving nominal significance (IVW P < 0.05) but not surviving FDR correction were classified as exploratory signals requiring independent replication and are discussed as such throughout this manuscript. When Cochran’s Q test indicated no significant heterogeneity (P > 0.05), we applied the fixed-effects IVW model. Conversely, when substantial heterogeneity was detected (Cochran’s Q test, P < 0.05), we employed the random-effects IVW framework. The weighted median, MR-Egger, and weighted mode approaches were implemented as sensitivity analyses to assess the robustness of IVW findings. Effect estimates are reported as odds ratios (ORs) with 95% confidence intervals. Where the primary IVW analysis yielded P < 0.05 but sensitivity methods were non-significant, we conservatively treated such associations as exploratory rather than confirmatory, given that consistent directional agreement across methods is required for stronger causal inference. These findings are explicitly flagged as requiring further validation. Scatter plots were generated to visualize MR results. This integrated analytical strategy enabled comprehensive assessment of potential causal pathways. Additionally, reverse MR analyses were performed to examine whether PsA liability causally influences immune cell profiles.

Sensitivity Analysis

To assess the reliability of our causal estimates, we implemented several sensitivity tests. Heterogeneity across instrumental variables was quantified using Cochran’s Q statistic, which evaluates consistency of SNP-specific causal effects. We examined potential violations of the exclusivity assumption by testing for horizontal pleiotropy through MR-Egger regression, wherein a non-zero intercept indicates that genetic variants influence the outcome through pathways independent of the exposure.²³ To confirm the directionality of causal inference and ensure that observed associations reflect effects of immune phenotypes on PsA rather than reverse causation or confounding, we applied the Steiger filtering method, which compares the variance explained in exposure versus outcome by each instrumental variable.²⁴ Instruments for which the variance explained in the outcome exceeded that in the exposure were excluded prior to final estimation, further guarding against spurious findings that could arise from the relaxed instrument selection threshold applied in this study.

Results

Exploratory Analysis of Immunocyte Phenotypes on PsA

The primary IVW analysis identified 26 immunophenotypes with a suggestive nominal association with PsA risk (P < 0.05, see Figure 2). Among these, 14 traits were nominally linked to heightened disease susceptibility (OR > 1), predominantly within B-cell and T-cell compartments. Potential phenotypes included: BAFF-R expression on IgD⁺CD24⁺ B cells; CCR7 on naïve CD4⁺ T cells; CD127 on CD45RA⁺ CD4⁺ T cells; BAFF-R on IgD⁻CD38⁺ B cells; IgD on IgD⁺CD38⁻ unswitched memory B cells; IgD on IgD⁺CD24⁺ B cells; and CD25 expression across multiple B-cell subsets (CD24⁺CD27⁺, memory, IgD⁺CD24⁺, unswitched memory, and IgD⁻CD38⁻ populations); CD39 on monocytes; along with total T-cell and lymphocyte counts. In contrast, 12 immune traits demonstrated inverse nominal associations with PsA risk (OR < 1, P < 0.05): side scatter area (SSC-A) parameters on CD4⁺ T cells and plasmacytoid dendritic cells; absolute cell counts of CD33dim HLA-DR⁺CD11b⁻ populations; quantitative and proportional measures of CD8dim natural killer T cells; CD8 expression on terminally differentiated CD8⁺ T cells; CD8dim NKT cell frequency within the T-cell pool; CD3 on resting CD4⁺ regulatory T cells; CD28⁺CD45RA⁻ and CD45RA⁺CD28⁻ CD8⁺ T-cell proportions; and absolute counts of both CD45RA⁻CD28⁻ and CD45RA⁺CD28⁻ CD8⁺ T-cell subsets. When subjected to multiple testing correction via the false discovery rate approach, none of these associations achieved statistical significance (all FDR > 0.05), indicating that these findings represent preliminary, hypothesis-generating signals.

Figure 2 Results from Forward MR Analysis Using Inverse Variance Weighting. This figure presents the outcomes of a forward MR analysis examining the association between various immune cell markers and the risk of psoriatic arthritis. The analysis was performed using the inverse variance weighting (IVW) method. Each row represents a specific immune cell type or marker (e.g., CD45RA–CD28–CD8+ T cell) with corresponding regression coefficients (b), odds ratios (OR), and 95% confidence intervals (CI). The nominal P-values indicate suggestive associations, with those below 0.05 defined as nominally significant (highlighted in bold). Importantly, these associations did not remain significant after False Discovery Rate (FDR) correction and should be interpreted as exploratory. The chart also includes error bars representing the 95% confidence intervals for each odds ratio. Data were derived from a large cohort study, and the analysis accounts for multiple immunological exposures to assess their role in the development of psoriatic arthritis.

Robustness Assessment of Forward Analysis

Sensitivity testing confirmed the stability of these nominally significant findings. For each of the 26 identified immune traits, Cochran’s Q-test revealed no substantial instrument heterogeneity (all P > 0.05), and MR-Egger intercept analyses detected no evidence of directional pleiotropy (all intercept P > 0.05), indicating that the IVW estimates were not meaningfully biased by detectable pleiotropic variants. Weighted median and MR-Egger point estimates were directionally consistent with the primary IVW results for the majority of associations, further supporting their internal robustness within this exploratory framework. Full sensitivity statistics are provided in Supplementary Files 1 and 2.

Bidirectional Relationship: PsA Effects on Immune Profiles

To examine whether disease liability potentially influences immune composition, we performed reverse MR analyses treating PsA as the exposure (Figure 3). Genetic predisposition to PsA showed preliminary positive nominal associations (OR > 1, P < 0.05) with five immune traits: CD45 expression on lymphocytes, total T cells, and HLA-DR⁺ CD4⁺ T cells; IgD⁺CD38dim B-cell frequency among B cells; and IgD⁺CD38dim B-cell proportion within total lymphocytes. Conversely, PsA demonstrated inverse nominal associations (OR < 1, P < 0.05) with 18 immunophenotypes. These included reduced HLA-DR expression on myeloid dendritic cells, CD33dim HLA-DR⁺CD11b⁻ populations, CD14⁺CD16⁺ monocytes, and HLA-DR⁺ CD4⁺ T cells. Additional traits showing negative associations comprised specific T-cell markers (CD8 on CD28⁺CD45RA⁺ CD8⁺ T cells, CD3 on CD28⁻ CD8⁺ T cells, CD4⁺CD8⁺ double-positive T-cell counts) and B-cell characteristics (CD25 on IgD⁺CD38dim B cells). Consistent with forward analyses, these associations remained nominally significant but did not survive FDR correction, indicating they represent preliminary signals requiring independent validation.

Figure 3 Results from Reverse MR Analysis Using Inverse Variance Weighting. Forest plot depicting odds ratios (OR) with 95% confidence intervals (CI) for the nominal associations between 23 immune cell phenotypes and disease susceptibility. Each row represents a distinct immune cell trait with corresponding effect estimates (b), number of single nucleotide polymorphisms (SNPs) used as instrumental variables, and P values. IVW, inverse variance weighting method. Immune cell traits with OR < 1 (left of vertical reference line) indicate protective factors, whilst those with OR > 1 (right of vertical reference line) indicate risk factors. Error bars represent 95% confidence intervals. All associations achieved nominal statistical significance (P < 0.05) but did not survive FDR correction.

Reverse Sensitivity Analyses

Robustness checks for the 23 reverse associations supported the validity of selected genetic instruments. Cochran’s Q-test indicated no significant heterogeneity across instruments (all P > 0.05), and MR-Egger intercept analyses found no systematic directional pleiotropy (all P > 0.05). Weighted median estimates were directionally concordant with IVW results for the majority of reverse associations, reinforcing the consistency of these preliminary signals. Detailed sensitivity statistics are presented in Supplementary Files 3 and 4.

Discussion

Our analysis applies a two-directional, two-sample Mendelian randomization approach to explore how immune cell characteristics nominally influence PsA susceptibility, and conversely, how disease predisposition potentially reshapes immune profiles across 731 distinct cellular features (Figure 4). Despite therapeutic advances centered on blocking the IL-23/IL-17 pathway, many patients continue to experience inadequate treatment response, pointing toward the putative involvement of immune mechanisms beyond this axis.^1,25 We uncovered preliminary putative evidence linking 14 immune traits to heightened PsA risk and 12 traits to reduced risk. These findings nominally shift attention from the predominant focus on T lymphocytes by suggesting that B-cell activation may play a previously underappreciated exploratory pathogenic role, while also aligning with the contribution of specific T-cell populations and demonstrating that PsA liability is nominally associated with altered myeloid cell function.

Figure 4 Schematic diagram of the conclusion. Abbreviations: IgD, immunoglobulin D; PsA, psoriatic arthritis.

Among our most striking observations is the suggestive association between activated B-cell populations and PsA liability. In particular, heightened BAFF receptor density on IgD⁺CD24⁺ B cells emerged as a notable candidate risk determinant. Although PsA is traditionally classified within the spondyloarthritis spectrum—disorders characterized primarily by innate and T-cell inflammation²⁶—our data align with a growing body of evidence implicating humoral immunity in disease chronicity. The putative significance of BAFF-R suggests that survival signals required to maintain potential auto-reactive B-cell clones could be critical upstream drivers. This offers a hypothetical rationale for the formation of ectopic lymphoid structures observed in the psoriatic synovium, which serve as local germinal centers for antibody production.^27–30 Importantly, although our results are similar to those of Jia et al (2025) on T cells,⁸ our work diverges and adds value in two fundamental ways. First, whereas Jia et al reported a lack of reverse causal effects, our bidirectional analysis identifies nominal evidence of PsA-driven immune remodeling—specifically the depletion of circulating pDCs and alterations in myeloid compartments—which better aligns with observed clinical phenomena of cellular recruitment to inflamed joints. Second, while prior work provided a broad-spectrum overview, we provide a deep-dive into the B-cell co-stimulatory axis (notably CD25 and BAFF-R signals). By framing these findings as a hypothesis-generating map rather than a definitive list, we offer a more nuanced interpretation of how genetic liability to B-cell activation might precede clinical diagnosis, a conceptual layer that was not the primary focus of the Jia et al study. Furthermore, recent molecular profiling has identified specific PsA endotypes characterized by a “B-cell rich” synovial signature, suggesting that B-cell depletion or BAFF blockade, strategies which have yielded equivocal results in unstratified trials, may warrant re-evaluation in these specific patient subgroups.^30,31

Our findings nominally highlight the potential pathogenic contribution of specific adaptive immune cell populations to PsA susceptibility. Suggestive evidence linked increased CCR7 expression on naïve CD4⁺ T cells and elevated CD127 on CD45RA⁺ CD4⁺ T cells with heightened disease risk. These genetic signals results align with transcriptomic observations by Ezeonyeji et al, who showed that phenotypically naïve T cells from psoriatic patients exhibit functional abnormalities, secreting disproportionate amounts of IFNγ and IL-22 that in turn stimulate keratinocyte production of the chemokine CXCL9.^32,33 Recent single-cell sequencing has further confirmed the clonal expansion of specific pro-inflammatory T-cell subsets in the PsA joint.³⁴ However, it is vital to distinguish these genetically proxied liabilities from transient cellular dynamics observed in functional studies. Furthermore, the nominal link between increased T-cell subsets and PsA consistent with the efficacy of existing T-cell-modulating therapies. For instance, Alefacept, which targets CD2 and reduces T-cell counts, has been shown to significantly improve joint symptoms and reduce synovial T-cell infiltration.³⁵ Similarly, Fludarabine treatment in PsA patients induces moderate synovial lymphopenia which correlates with clinical response, further cementing the pathogenic role of T-cell burden.^36,37 However, our identification of novel markers like CD127 (IL-7R) suggests that more refined targeting of IL-7 signaling pathways could offer superior specificity compared to broad lymphocyte depletion.

Conversely, we observed that 12 phenotypes, including SSC-A on plasmacytoid dendritic cells and absolute counts of CD33^dimHLA-DR⁺CD11b⁻ cells, were nominally linked to a reduced likelihood of PsA. Rather than implying that these cells are inherently “protective,” this negative correlation likely reflects the “drainage” or “tissue homing” hypothesis. As separate clinical studies have firmly established that circulating pDCs are reduced in the peripheral blood of PsA patients because they are actively recruited to the inflamed joint via chemotactic gradients, such as those mediated by CXCR3 and its ligands (CXCL10/CXCL9).^38–40 Indeed, pDCs are significantly enriched in PsA synovial fluid compared to osteoarthritis or healthy controls.^41–44 Our MR findings provide a complementary genetic proxy for this migration capacity, though they do not directly measure cellular flux. Therefore, the genetically predicted “lower” levels in blood may serve as a preliminary proxy for heightened capacity to migrate into tissue, mirroring the active consumption of these cells in the joint microenvironment.

However, several limitations must be acknowledged. First, the exposure GWAS relied on a relatively modest sample size (n = 3,757), which limited statistical power and necessitated the use of a relaxed P-value threshold for instrument selection. Although all instrumental variables (IVs) maintained F-statistics > 10, a degree of potential weak instrument bias cannot be entirely ruled out, which may influence the precision of our causal estimates. Second, due to the high-dimensional nature of the 731 immunophenotypes analyzed, the associations did not retain statistical significance after FDR correction. Therefore, these findings should be interpreted as exploratory and suggestive nominal associations requiring stringent validation in larger, independent cohorts. Third, the use of summary-level data inherently precludes the exploration of non-linear associations or stratified analyses (eg., by disease severity or medication history), which might mask endotype-specific immune drivers. Fourth, the genetic instruments employed in this study reflect constitutive, genetically determined expression levels rather than transient or context-dependent immune activation states; accordingly, caution is warranted when interpreting these findings as direct evidence of disease-relevant cellular dynamics, particularly during acute disease flares. Additionally, while we identify specific surface markers (eg., CD25, BAFF-R), the precise downstream effector functions of these subsets in the joint microenvironment remain to be fully characterized via functional assays. Importantly, the immune phenotypes identified here represent genetically proxied liability to altered surface marker expression or cell proportions, and should not be equated with functionally activated or pathologically expanded populations in the joint microenvironment. Any therapeutic inference drawn from these associations remains speculative and requires corroboration from functional and clinical studies.

Conclusion

In conclusion, our MR-based genetic atlas of immune cell phenotypes in PsA provides exploratory evidence for the involvement of the B-cell axis alongside traditional T-cell pathways. By identifying nominally associated B-cell traits (CD25⁺, BAFF-R⁺) as candidate contributors to disease risk, we have extended the presumed pathogenic situation beyond the T-cell/cytokine axis. It is important to emphasize that while these findings provide a genetically informed map of immune traits, all reported associations were nominal (P < 0.05) and did not survive FDR correction. Consequently, the identified pathways, including targeting B-cell co-stimulation and restoring regulatory cell populations, represent hypothesis-generating avenues rather than confirmed therapeutic strategies. Rigorous validation in larger cohorts and functional assays is necessitated to transition these preliminary signals into established clinical targets.