Abstract
Background:
Gut microbiota is involved in the development of psoriatic arthritis (PsA), but until now, there has been a lack of understanding of the PsA host–bacteria interaction.
Objectives:
To reveal the labels of gut microbiota in PsA patients and the species and functions related to disease activity.
Design:
Observational research (cross-sectional) with an exploratory nature.
Methods:
Metagenomics sequencing was used to analyze stool samples from 20 treatment-naïve PsA patients and 10 age-matched healthy individuals. All samples were qualified for subsequent analysis.
Results:
Compared with the healthy group, α-diversity was reduced in the PsA group, and β-diversity could distinguish the two groups. Two bacteria with high abundance and correlation with PsA disease activity were identified, Bacteroides sp. 3_1_19 and Blautia AF 14-40. In different functions, K07114 (calcium-activated chloride channel (CaCC) homolog) showed a positive correlation with PsA disease activity (disease activity in psoriatic arthritis, DAPSA) and Tet32 (an antibiotic-resistant gene), and carbohydrate-binding module family 50 was negatively correlated with erythrocyte sedimentation rate. A bacterial co-expression network associated with DAPSA was constructed. The network was centered on the bacteria in the Bacteroides genus, which formed a closely related network and were positively correlated with DAPSA. As another core of the network, K07114 was closely related to multiple bacteria in the Bacteroides genus and is also positively correlated with disease activity.
Conclusion:
The network composed of Bacteroides is associated with PsA disease activity, and its therapeutic value needs to be further explored. CaCCs may be a key channel for the interaction between Bacteroides and PsA-host.
Background
Psoriatic arthritis (PsA) is a psoriatic-related chronic inflammatory musculoskeletal disorder belonging to spondyloarthritides. 1 The development of PsA is involved in the interactions among genetic background, environmental factors, and the immune system. 2 The intestine is one of the main organs in the communication between the external environment and the human body, and keeping a delicate balance between the gut microbiota and the human immune system is crucial. It has been suggested that destabilization of balance may trigger the development of diseases.3,4
The immune inflammatory response-driven differentiation of T helper 17 (Th17) cells is a key point in the initiation and amplification of inflammation in PsA. 5 Gut bacteria modulate Th17 cells by producing bile acid metabolites. 6 Short-chain fatty acids (SCFAs) produced by gut bacteria fermenting carbohydrates have roles in processes such as regulating immune activation-related genes, the differentiation of T cells, and cytokine production. Moreover, SCFA pentanoate can effectively inhibit the proliferation of Th17 lymphocytes and the production of interleukin-17A (IL-17A) and can upregulate the expression of IL-10, which has the potential value in the treatment of inflammatory diseases. 7 HLA-B27 is one of the genes associated with PsA. Gut microbiota and their metabolites play a role in the development of HLA-B27-related diseases. 8
The reduced diversity of the gut microbiota has been found in PsA patients, with significant reductions in Coprococcus spp., Akkermansia, Ruminococcus, and Pseudobutyrivibrio. 9 But so far studies focusing on the gut microbiota of PsA patients are scarce and have been conducted at the genus level (by 16S rDNA sequencing) only. 10 We all know that a genus contains many different species that may play various roles in disease initiation and progression.11,12 Therefore, a more precise sequencing approach is required to deeply mine the gut microbiota composition into the species level in PsA patients.
Shotgun metagenomics is suitable for more accurate analysis of the gut microbiota. 12 Shotgun metagenomics, designed for complete DNA sequencing in specific environments, offers advantages such as wide coverage, accurate taxonomy (species level), and reliable function. This study aimed to explore the signatures of the gut microbiota of PsA patients at a greater sequencing depth and further analyze the correlation between these signatures and disease activity.
Methods
Subjects and sample collection
In this cross-sectional study, 20 treatment-naïve PsA patients and 10 age-matched healthy controls were enrolled. All the included PsA patients were aged 18–60 years and fulfilled the 2006 Classification criteria for Psoriatic Arthritis. The exclusion criteria were as follows: (1) use of antibiotics or probiotics in 1 month before collection of the stool sample; (2) long-term use of proton pump inhibitors, statins, metformin, and antipsychotics; (3) complicated with hepatitis B, tuberculosis, or other chronic infections; and (4) pregnancy or lactation. Healthy controls were recruited from people who came to our center for physical examination and met all the following: (1) No previous chronic disease history, (2) no chronic disease was found after physical examination, (3) no history of acute infection and antibiotic use in the past 1 month, and (4) no diarrhea in the past 1 month.
The stool containers were prepared, hands were thoroughly washed, and gloves were worn when collecting fresh stool specimens. The mid-piece of the sample was intercepted with a sterile stool sampler. One spoonful of fecal specimen was taken per tube with sterilized centrifuge tubes, and multitube backup was taken for each specimen. Fecal samples were kept in a −80°C freezer within 1 h after collection.
Metagenomic sequencing
Shotgun metagenomic sequencing based on high-throughput sequencing technology was utilized to obtain the complete images of microbial communities of all involved participants. DNA was extracted from the samples for library construction and library detection, and then sequencing based on the qualified libraries was performed with Illumina PE150 platform (Illumina, San Diego, CA, USA).
The sequencing data were used for post-information analysis. Quality control and host filtering of the raw data, as well as metagenome assembly and gene prediction, were performed to obtain a gene catalog. Functional annotation and abundance analysis of metabolic pathways (KEGG), carbohydrate enzymes (CAZy), and resistance genes (CARD) were conducted from the gene catalog.
Collection of clinical data of PsA patients
We collected demographics and clinical data including family history of psoriasis, symptoms, physical examinations, and laboratory tests. Disease history included hypertension, diabetes, major adverse cardiovascular and cerebrovascular events, daily medical treatments, and antibiotics (1–12 months). Laboratory tests included C-reactive protein (CRP), erythrocyte sedimentation rate (ESR), immunoglobulin G (IgG), immunoglobulin A (IgA), and immunoglobulin (IgM).
Scores assessing the disease activity of PsA, such as the Composite Psoriasis Disease Activity Index (CPDAI), disease activity in psoriatic arthritis (DAPSA), body surface area (BSA), and Psoriasis Area and Severity Index (PASI), were calculated. Central obesity was defined as a waist circumference (cm) greater than 85 cm in women and 90 cm in men.
The reporting of this study conforms to the Strengthening the Reporting of Observational Studies in Epidemiology statement. 13
Statistical analysis
R version 4.3.1 (R Foundation for Statistical Computing, Vienna, Austria) was used for statistical analysis, and packages such as “vegan,” “ggplot2,” “reshape2,” “weighted gene co-expression network analysis” (WGCNA), and “corrplot” were used for analysis of the data and mapping. The α-diversity, a comprehensive indicator of the number of microorganisms in a group, was calculated according to the Shannon‒Wiener index and differences between the two groups were compared by Student’s t-test, with a significant difference at p < 0.05. The β-diversity, a measure of the differences in the number and distribution of each species between groups, was displayed with principal coordinate analysis (PCoA), which was calculated based on the Bray‒Curtis distance, and the permutational multivariate analysis of variance (PERMANOVA) method was used to compare differences between two groups (permutations = 999). The top 10 bacteria by relative abundance at the family, genus, and species levels were demonstrated using stacked plots. The top 20 functions (level 2 of KEGG) were demonstrated using stacked plots. The correlation of the abundances of the top 10 antibiotic-resistant genes (ARGs) with individual samples was visualized using Circos.
Metastats analysis was used to compare taxonomic differences between the two groups at the species level according to the relative abundance of each sample, with a significant difference at q < 0.05 (after correction for p by the Benjamini and Hochberg method). Linear discriminant analysis (LDA) effect size (LEfSe) analysis was used to compare functional differences based on annotated genes in the CAZy, KEGG, and CARD databases between the two groups, with a significant difference at p < 0.05 and LDA score > 2.
The Spearman method was used to calculate the correlation between two variables, and the correlation coefficient r represents the correlation size and direction, and p < 0.05 was considered statistically different (Benjamini and Hochberg method correction).
WGCNA algorithm, a powerful systems biology approach, was used to extract bacterial co-expression modules from all differential bacteria (species level, PsA vs control). The co-expression modules were constructed by the WGCNA package in R. Cytoscape_v3.9.1 (Shannon et al., 2003) was used to draw network maps of key bacteria, clinical indicators, and differential functions.
Results
Demographics of PsA patients and healthy controls
This study included 20 PsA patients and 10 age-matched healthy controls. For PsA patients, the courses of psoriasis and arthritis were, respectively, 6.4 (0.0–19.7) years and 11.0 (3.8–24) months. There were no significant differences in age, sex, body mass index (BMI), waist circumference, smoking, or alcohol consumption between PsA patients and healthy controls, although PsA patients had a greater family history of psoriasis (40.0% vs 0.0%) (Table 1).
Demographic and clinical data of PsA patients and healthy Ctr.
MACCE is defined as one or more of the following conditions: acute myocardial infarction, cardiogenic shock or acute heart failure, unstable angina, and stroke. In daily medical treatments, one patient took nifedipine for a long time. Antibiotics (1–12 months), antibiotics were administered 1–12 months before inclusion in this study, PsA (five cases, β-lactam antibiotic; one case, quinolones antibiotic), Ctr (two cases, β-lactam antibiotic), all medication time is less than or equal to 3 days.
BSA, Body Surface Area; CPDAI, Composite Psoriasis Disease Activity Index; CRP, C-reactive protein; Ctr, controls; DAPSA, disease activity in psoriatic arthritis; ESR, erythrocyte sedimentation rate; MACCE, major adverse cardiovascular and cerebrovascular event; NA, not applicable; PASI, Psoriasis Area and Severity Index; PsA, psoriatic arthritis.
Metagenomic sequencing revealed a taxonomic signature of the gut microbiota in PsA patients
Microbial α-diversity calculated by the Shannon index at the species level was significantly lower in PsA patients than in healthy controls (1.55 ± 0.21 vs 1.73 ± 0.21, p = 0.038) (Figure 1(a)). PCoA based on the Bray–Curtis distance (β-diversity) showed bacteria separation of PsA patients from controls (PERMANOVA, p = 0.039, permutations = 999) (Figure 1(b)). The top 10 bacteria in terms of relative abundance in two groups are displayed in stacked plots at the family, genus, and species levels (Figure 1(c)–(e)). Metastats analysis revealed 1254 differential bacteria at the species level between two groups, with the top 10 differential species based on abundance displayed in boxplots (Blautia sp. AF14-40, Bacteroides sp. 3_1_19, Blautia sp. SG-772, Bacteroides sp. 2_2_4, Desulfotomaculum sp. OF05-3, Thermophilibacter provencensis, Lachnospiraceae bacterium AM40-2BH, Bacteroides sp. 3_1_23, Bacteroides sp. HF-5092, and Bifidobacteriaceae bacterium MCC01984) (Figure 1(f)).
Taxonomic differences of gut microbiota between PsA patients and controls. (a) α-Diversity: according to the Shannon–Wiener index. (b) β-Diversity: PCoA was calculated based on the Bray–Curtis distance, and the PERMANOVA method was used to compare differences. (c) Top 10 taxonomic communities at family. (d) Top 10 taxonomic communities at genus. (e) Top 10 taxonomic communities at species levels. (f) Species differences between PsA patients and controls: metastats analysis was used, with a significant difference at q < 0.05.
Metagenomic sequencing revealed a functional signature of the gut microbiota in PsA patients
We used the CAZy, KEGG, and CARD databases for functional annotation to explore the microbial function in PsA patients. The PCoA based on gene abundance at level 2 of the CAZy database showed functional community separation of PsA patients from controls (PERMANOVA, p = 0.046, permutations = 999) (Figure 2(a)). These pathways were annotated in the KEGG database at level 2, and the top 20 ranked by relative abundance are displayed in a stacked plot (Figure 2(b)). A total of 193 ARGs were annotated using the CARD database, and the top 10 ARG relative abundances and their corresponding fecal samples are shown in Figure 2(c).
Functional differences in gut microbiota between PsA patients and controls. (a) β-Diversity: PCoA was calculated based on the Bray–Curtis distance, and the PERMANOVA method was used to compare differences (CAZy, level 2). (b) Top 20 pathways ranked by relative abundance in the KEGG database (level 2). (c) Top 10 ARO relative abundance and their corresponding fecal samples (CARD). LEfSe analysis was used to compare functional differences, with a significant difference at p < 0.05 and LDA score > 2. (d) LEfSe analysis in CAZy. (e) LEfSe analysis in KEGG. (f) LEfSe analysis in CARD.
LEfSe analysis was used to explore the functional signature of the gut microbiota in PsA patients. In the CAZy analysis, glycoside hydrolase family 43 (GH43) was upregulated, while carbohydrate esterase family 4 (CE4) and carbohydrate-binding module family 50 (CBM50) were downregulated in PsA patients (Figure 2(d)). In the KEGG analysis, K02004 (putative ABC transport system permease protein), K01190 (beta-galactosidase), K05349 (beta-glucosidase), K01897 (long-chain acyl-CoA synthetase), K01187 (alpha-glucosidase), and K07114 (calcium-activated chloride channel (CaCC) homolog) were upregulated, while K06131 (cardiolipin synthase A/B), K07729 (putative transcriptional regulator), and K01005 (polyisoprenyl-teichoic acid-peptidoglyteichoic acid transferase) were downregulated in PsA patients (Figure 2(e)). In the CARD analysis, tet32 and Ccol_ACT_CHL were upregulated, while vanWG, vanYA, KpnEF, emrR, C-reactive protein (CRP), marA, and ANT6_la were downregulated in PsA patients (Figure 2(f)).
Bacteria and functions associated with PsA disease activity
The correlations between the differential bacteria (top 10) and differential functions (KEGG, CARD, CAZy) and PsA disease activity were shown to explore the macro–micro relationship (Figure 3(a)–(d)). The relative abundance of Bacteroides sp. 3_1_19 showed a significantly positive correlation with PsA disease activity indicators, including CPDAI, DAPSA, and ESR, as well as with serum IgA (Figure 3(a)). The relative abundance of Blautia sp. AF14-40 was significantly negatively correlated with DAPSA, PASI, and ESR in PsA patients (Figure 3(a)). K07114 (CaCC homolog) showed a significantly positive correlation with PsA disease activity indicators, including DAPSA and ESR (Figure 3(b)). Tet32, an ARG, was negatively correlated with ESR (Figure 3(c)). CBM50 was negatively correlated with ESR (Figure 3(d)).
Correlation between differential bacteria (top 10 species), differential functions (KEGG, CARD, CAZy), and clinical indicators. (a) Differential bacteria (top 10 species) and clinical indicators. (b) Differential functions (KEGG) and clinical indicators. (c) Differential functions (CARD) and clinical indicators. (d) Differential functions (CAZy) and clinical indicators. The area of the circle represents the absolute value size of the correlation coefficient (r) value. The color of the circle represents the positive/negative correlation and the depth of the color represents the size of the correlation coefficient (r) value. p Values were corrected by the FDR method in multiple testing.
Bacterial co-expression network associated with PsA disease activity
Based on all differential bacteria, four modules (“MEturquoise,” “MEblue,” “MEbrown,” “MEgrey”) were established by the WGCNA algorithm. The heatmap of the correlation between WGCNA modules and clinical indicators (Figure 4(a)) showed that “MEbrown” was the most strongly correlated with DAPSA. The bacterial co-expression network was constructed based on “MEbrown” and DAPSA score (Figure 4(b)). The network is centered on the bacteria in the Bacteroides genus. Several bacteria in the Bacteroides genus form a closely related network, which is significantly positively correlated with the DAPSA (Figure 4(b)). K07114 (CaCC homolog) was also included in the co-expression network because the above analysis found that it was significantly associated with DAPSA. K07114 was another core of this network, which was closely related to multiple bacteria of the Bacteroides genus and positively correlated with DAPSA.
Bacterial co-expression networks associated with PsA disease activity. (a) Correlation of each co-expression module (WGCNA) with clinical information. The numbers in the box represent Spearman’s correlation coefficient (p-value). (b) Co-expression network was constructed according to the bacteria associated with both the WGCNA MEbrown module and the DAPSA score. DAPSA-related differential functions (K07114, Ca-activated chloride channel homolog) were added to the network. Circles represent bacteria (species level), rectangles represent DAPSA, and hexagons represent microbial functions. The circle size represents the relative abundance level of bacteria, and the filled color of the circle represents the phylum to which the bacteria belong. Blue rings represent the Bacteroides genus. The line between two points represents that there is a significant correlation between two points, red line is a positive correlation, blue line is a negative correlation.
Discussion
Shotgun metagenomics enables the studies of the gut microbiota not only down to the species level but also more reliable functional annotation. 12 Microbes, a second genome in humans, play important roles in immune regulation, metabolism, nutrient intake, endocrine regulation, and numerous diseases. The gut, which represents the body’s largest immune organ, stores a diverse array of immune cells in its lamina propria and lymphoid tissues. The gut is also a major site of mutual training between the human immune system and the microbiota. 14 Gut microbiota has been reported to be associated with the development of PsA, but all these studies focused only on the genus level instead of the species level. Moreover, to the best of our knowledge, no general conclusions on bacteria–host interactions can be drawn at the genus level.10,15 Shotgun metagenomics can provide a more precise characterization of complete DNA by sequencing all DNA within a certain environment.
Shotgun metagenomics was used in this study to explore taxonomic and functional differences in the gut microbiota between PsA patients and healthy controls. We showed a reduced diversity of the gut microbiota in patients with PsA, which is consistent with previous reports. 9 The previous data regarding the Bacteroides abundance in psoriatic patients were inconsistent, being elevated in some studies, while decreased in other studies. 10 These conflicting results may be related to insufficient depth of 16S rDNA sequencing. In the present study, we further revealed several significantly elevated species in the Bacteroides genus (Bacteroides sp. 3_1_19, Bacteroides sp. 2_2_4, Bacteroides sp. 3_1_23, and Bacteroides sp. HF-5092) in PsA patients compared with healthy controls. It is worth noting that the abundances of Bacteroides sp. 3_1_19 and Bacteroides sp. 2_2_4 were not only higher but also significantly positively correlated with DAPSA. Although Shapiro et al. reported the elevated abundance of Blautia in psoriasis patients, 16 we further revealed opposite changes at the species level in this genus (Blautia sp. AF14-40 elevated and Blautia sp. SG-772 decreased) in PsA patients. Different bacteria of the Blautia genus may play completely opposite roles in human health. 15 Our results are consistent with this view and provide some evidence. Our samples are from feces. If these differential bacteria have the potential to become therapeutic targets in the future, we suggest that these bacteria are located in the colon of the gastrointestinal tract.
Functional annotation by shotgun metagenomics has emerged as a robust method to identify aggregate functions encoded by community constituent microbes. 17 These homologous genes are involved in signaling pathways, including carbohydrate metabolism, lipid metabolism, as well as glycan biosynthesis and metabolism. Similar to our findings, Chen et al. 18 found that metabolic pathway genes involved in bacterial chemotaxis and carbohydrate transport were overrepresented in PsA patients. SCFAs serve as key points for host–bacteria interactions and can regulate inflammatory cytokines and T cells. 19 Recently, it was found that the gut microbiota can be restored in ankylosing spondylitis patients treated with tumor necrosis factor-alpha (TNF-α) inhibitors, and this restoration was achieved most likely by affecting SCFA-producing bacteria. 20 More specifically, our results showed enrichment of the K01897 ortholog. K01897 ortholog encodes a long-chain acyl-CoA synthetase that plays an important role in fatty acid biosynthesis and degradation. Carbohydrate metabolism is extremely closely linked to fatty acid metabolism. 21 In this study, we found that the abundance of two genes encoding enzymes acting on carbohydrates was decreased in PsA patients, one encoding CE4 and another encoding CBM50. ARG contamination has become a major public health problem worldwide. We also explored the ARGs that may be associated with PsA. ARGs are associated with the progression of diabetes and are negatively correlated with butyrate, a metabolite of the gut bacteria. 22 At present, the association of ARGs with immune diseases is far from being thoroughly studied. For the first time, we showed that two ARGs were enriched in patients with PsA, but further studies are warranted to explain their significance in PsA patients.
We found two bacteria with high abundance and correlated with the disease activity of PsA, Bacteroides sp. 3_1_19 and Blautia AF 14-40. We also found a function called K07114 (CaCC homolog), which was upregulated in PsA patients and positively correlated with disease activity. At present, we cannot fully understand the biological significance of these new findings.
However, it was interesting that Bacteroides sp. 3_1_19, Blautia AF 14-40, and K07714 were in the same co-expression network. The network was centered on the bacteria in the Bacteroides genus. Several bacteria in the Bacteroides genus formed a closely related network, which was significantly positively correlated with the DAPSA. Bacteroides genus, as a dominant component of the Bacteroidetes phylum, most studies support its beneficial effects against obesity, metabolic disorders, and inflammation. 23 Several bacteria in the Bacteroides genus have been identified with the ability to relieve lipopolysaccharide (LPS)-induced IL-8 release, a model for evaluating inflammation. 24 The co-expression network in our results showed that several members of the Bacteroides genus were closely formed a bacterial community and were positively associated with inflammation levels in PsA patients. Whether this network has therapeutic value needs to be confirmed by randomized controlled trials.
K07114 was closely related to multiple bacteria of the Bacteroides genus and was one of the cores of this network. This result suggested that the Bacteroides may affect PsA disease activity through CaCCs. The gut epithelial barrier and vascular barrier are the main components of the gut barrier and are the key barriers that prevent microorganisms and their metabolites from reaching the systemic circulation. 25 Activation of CaCCs under inflammatory conditions aggravates epithelial barrier dysfunction in intestinal epithelial cells and vascular permeability, promoting microbiota-driven gut barrier leakage.25,26 Proinflammatory compounds produced by bacteria, such as LPS, are transferred from the “leaky gut” into the blood, leading to systemic inflammation, which then reaches the joints, with consequences such as osteoarthritis, rheumatoid arthritis, and spondyloarthritis (including PsA). 27
There are some limitations in this article. First, the sample size is small. The included subjects were treatment-naïve PsA patients, which could not reflect the potential impact of treatment on gut microbiota. Second, this study adopted a cross-sectional design and exploratory method, with no predefined specific hypothesis and no predefined statistical analysis plan. Therefore, no causal relation between the microbial findings and PsA disease can be determined in this study. These findings are also not representative of general PsA patients. Third, PsA patients have an increased risk of comorbidities, and it is difficult to completely rule out the influence of comorbidities. Overall, the findings of this study need to be confirmed in a greater prospective cohort.
In the future, there are still some problems worth exploring in the interaction between PsA and gut microbiota. It has been highlighted that certain treatments in spondyloarthritis, such as TNF inhibitors, have led to the expansion of specific bacterial species. 28 The impact of new treatments for PsA, such as IL-17 inhibitors or Janus kinases (JAK) inhibitors, on the gut microbiota of PsA deserves further study. In addition, diet and probiotic intake appear to play a role in the microbiota in spondyloarthritis (SpA). 29 It will also be of practical significance to include the analysis of dietary factors in future studies.
Conclusion
The gut microbiota of PsA had its signature, which was manifested by significant changes in some species (such as Bacteroides sp. 3_1_19 and Blautia AF 14-40), and some functions were up/downregulated (such as K07114, Tet32, and CBM50). It was noteworthy that these changes were associated with the disease activity of PsA. Bacteroides sp. 3_1_19, Blautia AF 14-40, and K07114 existed in the same co-expression network. The network composed of Bacteroides is associated with PsA disease activity, and its therapeutic value needs to be further explored. CaCCs may be a key channel for the interaction between Bacteroides and PsA-host.
