Abstract
Objective. To assess whether gene expression signatures associated with rheumatoid arthritis (RA) before pregnancy differ between women who improve or worsen during pregnancy, and to determine whether these expression signatures are altered during pregnancy when RA improves or worsens.
Methods. Clinical data and blood samples were collected before pregnancy (T0) and at the third trimester (T3) from 11 women with RA and 5 healthy women. RA disease activity was assessed using the Clinical Disease Activity Index (CDAI). At each timepoint, RA-associated gene expression signatures were identified using differential expression analysis of RNA sequencing profiles between women with RA and healthy women.
Results. Of the women with RA, 6 improved by T3 (RAimproved), 3 worsened (RAworsened),and 2 were excluded. At T0, mean CDAI scores were similar in both groups (RAimproved 11.2 ± 9.8; RAworsened 13.8 ± 6.7; Wilcoxon rank-sum test: P = 0.6). In the RAimproved group, 89 genes were differentially expressed at T0 (q < 0.05 and fold change ≥ 2) compared to healthy women. When RA improved at T3, 65 of 89 (73%) of these genes no longer displayed RA-associated expression. In the RAworsened group, a largely different RA gene expression signature (429 genes) was identified at T0. When RA disease activity worsened at T3, 207 of 429 (48%) genes lost their differential expression, while an additional 151 genes became newly differentially expressed.
Conclusion. In our pilot dataset, pre-pregnancy RA expression signatures differed between women who subsequently improved or worsened during pregnancy, suggesting that inherent genomic differences may influence how pregnancy affects disease activity. Further, these RA signatures were altered during pregnancy as disease activity changed.
Rheumatoid arthritis (RA) is a systemic inflammatory disease that leads to significant disability resulting from pain and swelling in inflamed joints and from joint destruction. To date, there is no cure. Pregnancy is known to have disease-modifying properties1,2,3,4 on RA, with a significant proportion of women experiencing an improvement in disease activity during pregnancy while in others, the disease may remain unchanged or even worsen. Even though there are medications, including some traditional and biologic disease-modifying antirheumatic drugs that are considered safe for use in pregnancy,5,6 many women with RA prefer to stop taking medications during pregnancy. However, because there are no known biomarkers at present to predict who is likely to improve or worsen, or whose RA will remain unchanged during pregnancy, these women are hesitant to plan a pregnancy because they do not know whether their disease will worsen if they stop taking medications in order to try to conceive.
Several case-control studies based on gene expression data from microarrays7,8,9,10,11,12 or RNA sequencing (RNA-seq) have been conducted to investigate gene expression signatures associated with RA.13 However, gene expression studies that have been conducted in the context of RA pregnancy have not examined RA-associated expression signatures in the nonpregnant state due to pre-pregnancy samples not being available.14,15,16,17 It is thus unknown whether the pre-pregnancy RA expression signature can be used to predict whether RA will subsequently improve or worsen during pregnancy. Further, given that gene expression is a dynamic process, it is possible that the RA-associated gene expression signature may be altered during pregnancy. Genes modulating disease activity during pregnancy may show altered expression when disease activity changes over time, as their expression may either no longer be associated with RA when disease activity is low or in remission during pregnancy, or additional genes may show RA-associated expression when RA worsens during pregnancy. However, the influence of pregnancy on the RA gene expression signature, if any, has not been investigated.
In our present study, we have used our unique prospective pilot pregnancy cohort of women with RA and healthy women that includes a pre-pregnancy timepoint18,19 as a case-control dataset to examine gene expression signatures associated with RA at the pre-pregnancy baseline. We hypothesized that the baseline RA-associated gene expression signature among women who subsequently improved during pregnancy differs from that of women who worsened during pregnancy. We also evaluated a second hypothesis that the gene expression signature associated with RA at pre-pregnancy is altered during pregnancy when disease activity improves or worsens.
METHODS
Study subjects. Healthy women and women with RA of Danish descent who were planning a pregnancy were recruited, enrolled in our pregnancy cohort in Denmark, and prospectively followed, as previously described.18 A subset of 11 women with RA and 5 healthy women from this cohort, on whom we reported longitudinal changes in expression,19 was included in the present study. The women with RA fulfilled the 1987 revised American College of Rheumatology criteria for the disease.20 The study was approved by the Ethics Committee for Region Hovedstaden (Denmark), the Danish Data Protection Agency, and the Children’s Hospital Oakland Research Institute institutional review board (IRB number: 2009-073). All subjects provided written informed consent prior to enrollment. Data and samples were collected as previously described.19
Assessment of RA disease activity. RA disease activity was assessed using the Clinical Disease Activity Index (CDAI)21 because it does not include acute-phase reactants such as C-reactive protein (CRP), whose levels are known to fluctuate during pregnancy;22,23 further, acute-phase reactants do not contribute much information on top of what is provided by the CDAI.24 The change in CDAI (ΔCDAI) from before pregnancy (T0) to the third trimester (T3) was used to determine whether disease activity improved or worsened. Patients were categorized as having improved by T3 (RAimproved) if their ΔCDAI fit the criteria for a minimum clinically important difference (MCID) based on baseline (T0) disease activity; ΔCDAI values of 12, 6, and 1 were used as thresholds when disease activity at T0 was high, moderate, or low, respectively.25 Those women with an increase in CDAI from T0 to T3 satisfying the MCID criteria for worsening of disease activity were included in the “worsened” subset, referred to as RAworsened.
RNA sequencing and bioinformatic analyses. RNA extractions, processing, and sequencing were performed as originally described.18 Pseudoalignment of the demultiplexed raw sequence reads (FASTQ format) to the Ensembl reference human GRCh38 transcriptome assembly (release 98) and quantification of transcript abundances were performed using kallisto (v0.43.0).26 BioMart annotations were used to combine transcript-level counts into gene-level estimates; counts were summed by Ensembl gene IDs. Gene IDs that mapped to patches or alternate haplotypes rather than to the primary reference sequence were excluded to avoid duplication. Pseudogenes, genes without annotations, and genes with very low read counts (< 1 transcript/million) in at least 25% of all samples were filtered out. Any globin and ribosomal RNA transcripts still present were also filtered out. To adjust for variable sequencing depths across samples, the gene-level counts were normalized using the trimmed mean of M values algorithm as implemented in the edgeR package (v3.26.8; http://bioconductor.org).27,28 To assess batch effects, normalized counts from pairs of technical replicates were plotted, and outliers were filtered out to achieve a Pearson correlation of at least 95% between replicates.
Case-control differential gene expression analysis. To identify gene expression signatures associated with RA at the T0 baseline, cross-sectional differential expression analysis was performed using edgeR (v3.26.8),27 comparing normalized T0 gene-level counts between each RA subset (RAimproved or RAworsened) and healthy women. In each analysis, a negative binomial distribution was used to handle the overdispersion in RNA-seq gene counts. Differential expression was tested using generalized linear model (GLM) likelihood ratio tests (LRT), and differences between women with RA and healthy women were assessed using the contrast argument of the glmLRT function in edgeR. Correction for multiple testing was performed using the false discovery rate method.29 A q-value threshold of 0.05, in combination with a fold change (FC) of at least 2, was used to assess significance. To determine whether the pre-pregnancy expression signature changed when RA improved or worsened during pregnancy, the differential expression analysis was repeated using data from the same women (RAimproved or RAworsened vs healthy women) at the T3 timepoint.
Functional analysis. Differentially expressed genes were analyzed for overrepresentation of Gene Ontology (GO) categories using a hypergeometric test implemented in the WEB-based GEne SeT AnaLysis Toolkit.30 A significance threshold of q < 0.05 was used to define enrichment. Cystoscape31 was used for functional annotations and visualization of protein interactions documented in the STRING database.32
RESULTS
Study subjects. Of the 11 women with RA, 6 improved by T3 while 3 worsened, based on MCID thresholds. Two women were excluded because even though their disease activity improved during pregnancy, one was already in remission at T0, and the ΔCDAI value for the other (ΔCDAI = 2.7) did not meet the MCID threshold of 6 for moderate baseline disease activity. The changes in disease activity scores from T0 to T3 were significantly correlated between the CDAI and the Disease Activity Score in 28 joints based on CRP 3 (DAS28-CRP3; Pearson correlation = 85%, P = 0.004). The average ages at conception were as follows: 28.9 ± 6.0 years for RAimproved, 33.2 ± 1.9 years for RAworsened, and 31.2 ± 5.7 years for the healthy women. The women who improved had a shorter disease duration than those who worsened (RAimproved [mean ± SD]: 6.5 ± 4.2 yrs; RAworsened: 8.9 ± 1.1 yrs), although this difference was not statistically significant. Medications taken by the women with RA at each timepoint are shown in Table 1. While mean disease activity (CDAI scores) at baseline did not differ significantly between the 2 RA subsets (RAimproved: 11.2 ± 9.8; RAworsened: 13.8 ± 6.7; Wilcoxon rank-sum test: P = 0.6), the mean values at T3 differed significantly (RAimproved: 2.2 ± 1.3; RAworsened: 31.7 ± 15.1; Wilcoxon rank-sum test: P = 0.02).
The RA gene expression signature at the T0 (pre-pregnancy) baseline. A total of 89 genes were differentially expressed (q < 0.05; FC ≥ 2) between the 6 RAimproved and 5 healthy women (Figure 1A; Supplementary Table 1, available with the online version of this article). The genes that were overexpressed (n = 44) in RA (e.g., C4BPA, CAMP, CD177, CRISP3, HLA-DQA2, MMP8, OLFM4, ORM1, S100A12) as well as those that were underexpressed (n = 45; e.g. CMPK2, HERC5, IFI44, IFI44L, IFITM3, IL1RL1, IL5RA, MX1, OAS1, OAS2, OAS3, SIGLEC1) were enriched in various immune-related GO biological processes, as shown in Table 2, and in Reactome pathways relating to interferon (IFN) signaling (q = 5.9 × 10−6), antiviral mechanism by IFN-stimulated genes (q = 7.0 × 10−5), and p130Cas linkage to mitogen-activated protein kinase signaling for integrins (q = 1.4 × 10−2), among others. A large proportion (52%) of the 89 genes differentially expressed between RAimproved and healthy women at the T0 baseline encode proteins that are functionally related, as shown by protein networks based on the STRING database32 in Cytoscape (Figure 2). Most of the underexpressed genes formed a tight cluster, distinct from the overexpressed genes.
A total of 429 genes were differentially expressed (FC ≥ 2; q < 0.05) between RAworsened and healthy women at T0 (Figure 3A). This gene expression signature largely differed from the one identified above (RAimproved vs healthy women). Only 19 of the 429 genes overlapped with those differentially expressed between RAimproved and healthy women at T0, with the majority demonstrating similar expression patterns in both RA subgroups compared to healthy women (overexpressed: OLFM4, UBB, ORM1, SEPTIN3, KRT1, TUBB2A; underexpressed: IL1RL1, IGLC3, IGLV2-14, PF4V1, FADS2, NKX3-1; Supplementary Table 2, available with the online version of this article). HLA-DRQA2, on the other hand, was 3-fold overexpressed among the RAimproved group and 3-fold underexpressed among the RAworsened women compared to the healthy women.
The RA gene expression signature is altered when RA improves or worsens during pregnancy. When disease activity improved by T3, most of the baseline RA signature genes identified among RAimproved women (65 of 89; 73%) were no longer differentially expressed between the RAimproved and healthy women (Figure 1B; Supplementary Table 1, available with the online version of this article). Of note, there were 24 genes that remained differentially expressed between RAimproved and healthy women at T3 (e.g., C4BPA, HLA-DQA2, IGLC3, IL5RA, MAOA, OLIG2, PTGDR2, TUBB2A), and an additional few (n = 27) became newly differentially expressed (e.g., ADORA3, CYP27A1, DSC1, RAP1GAP, SCL29A1).
When disease activity worsened by T3, 207 of the 429 genes (48%) differentially expressed at the T0 baseline lost their differential expression, while an additional 151 genes became newly differentially expressed (Figure 3B).
DISCUSSION
In the present study, our goal was to examine whether the pre-pregnancy RA gene expression signature differs between women who subsequently improved during pregnancy and those who worsened during pregnancy. We also examined whether the pre-pregnancy RA gene expression signature was altered in any way during pregnancy, when disease activity improved or worsened.
In our pilot dataset, even though mean disease activity was similar between the RAimproved and RAworsened groups at the pre-pregnancy baseline, there was very little overlap in the sets of genes showing RA-associated expression within each of the 2 groups. The few genes that overlapped between the 2 RA expression signatures included some that have previously been implicated in RA, such as IL1RL1,33 ORM1,34 KRT1,35 and HLA-DQA2.36 Although HLA-DQA2 expression was associated with RA in both subsets, this gene demonstrated contrasting expression patterns; it was 3-fold overexpressed in the RAimproved group and 3-fold underexpressed in the RAworsened group, both compared to the same set of healthy women. While increased HLA-DQA2 expression has been reported in RA,37 a negative correlation has also been found between expression levels in synovial tissue fibroblast cells and Health Assessment Questionnaire scores.38 Nonetheless, the significance of these contrasting expression patterns in the 2 groups of women with RA is not entirely clear. The RA expression signature identified among the RAimproved women included many additional genes whose expression and/or methylation patterns have previously been associated with RA, such as S100A12,39 CRISP3,40 MMP8,41 and CAMP.12 Of interest, the IFN-inducible genes IFI44, IFI44L, CMPK2, HERC5, MX1, SIGLEC1, OAS1, OAS2, and OAS3 were also part of the baseline RA expression signature among the RAimproved women. Compared to healthy women, these genes were underexpressed in this RA subset at baseline, as we previously reported,19 in contrast to other studies of RA and other autoimmune conditions.9,42,43 While our pre-pregnancy RA signatures also included many genes that did not overlap with RA signatures from previous case-control studies, results were inconsistent even across those previous studies. This could be attributed to a number of factors, including the source of RNA (whole blood [our study] vs peripheral blood mononuclear cells [PBMCs]),9,11,14,15,16,17 synovial fibroblasts,13 or neutrophils;44 differences in gene expression technology used (RNA-seq [our study] vs microarrays);9,11,14,15,16,17 and patient sample and sex ratio (only women, most of whom experienced improvement of RA disease activity during pregnancy [our study] vs women and men).9,11,13,44
In the present study, we also examined the influence of pregnancy on RA-associated gene expression signatures. This had not previously been reported since pre-pregnancy samples were not available in prior studies.14,15,16,17 We observed a dilution of the baseline RA signature during pregnancy, with the majority (73%) of signature genes showing similar expression in both women with RA and healthy women by the third trimester, when RA improved. These results are consistent with those of a previous study demonstrating minimal differences in PBMC expression profiles between women with RA and healthy women in the third trimester.15 While our study design did not allow us to determine whether the loss of association with RA in the third trimester is specific to pregnancy, it is plausible that as the expression profiles of the women with RA undergo pregnancy-induced changes, they start to resemble those of healthy women, as we observed, and these changes are accompanied by an improvement of the disease. On the other hand, when RA worsened during pregnancy, many genes (n = 151) demonstrated new expression patterns that became associated with RA, as would be expected when the “case” and “control” groups become phenotypically more different from each other.
In a previous study comparing the DAS28 based on erythrocyte sedimentation rate and DAS28-CRP, with and without patient global scores, during pregnancy, the DAS28-CRP3 was found to perform better.45 However, since even CRP levels have been shown to fluctuate during pregnancy,22,23 the DAS28-CRP3 does not represent a gold standard for use in pregnancy. Other measures of disease activity that do not include acute-phase reactants, such as the CDAI, have not been assessed for their performance during pregnancy. In the absence of a gold standard, we chose to use the CDAI to assess disease activity before and during pregnancy because acute-phase reactants have been shown to add little information on top of the clinical variables already included in the CDAI.24 Additionally, the CDAI is more stringent than the DAS28 when assessing improvement of disease activity; it has been shown that patients can satisfy DAS28 remission criteria while still having active disease.46
In our study, the availability of pre-pregnancy and pregnancy data from the same women who improved or worsened during pregnancy enabled us to investigate pre-pregnancy gene expression signatures between the 2 groups, as well as the effect of pregnancy on those expression signatures. The use of RNA-seq technology to assess gene expression was another strength. However, our study has some limitations. First, given that this is a pilot study, sample sizes were small. Nonetheless, patterns emerged that are supported by previous literature and thus further investigations in a larger sample are warranted. We did not examine proportions of cell types between disease states (RA vs healthy) and/or across timepoints because our goal was to identify overall systemic gene expression changes, resulting from altered expression of specific genes or from differences in cell proportions. Because we used total RNA from whole blood, expression profiles of neutrophils may have dominated the observed expression patterns. It is also possible that medications taken by the women with RA before pregnancy may have influenced the results. However, due to sample size limitations, we could not assess if this was the case, and we were also unable to adjust for dosage and/or specific medications.
In conclusion, we report novel findings that women with RA who improved during pregnancy demonstrated differences in pre-pregnancy RA-associated gene expression compared to women who worsened. These differences in pre-pregnancy RA expression signatures suggest that inherent genomic differences between women with RA may influence how pregnancy alters disease activity. Our findings that RA-associated gene expression signatures are altered during pregnancy when disease activity changes are also novel. Additional investigations in larger datasets are warranted to corroborate these preliminary findings and to identify novel drug targets and/or biomarkers of disease activity.
ACKNOWLEDGMENT
We are immensely grateful to the study subjects for their participation in the study. We acknowledge our gratitude toward Dr. Hanne Kjærgaard for all her efforts in setting up the logistics for this study in Denmark; Dr. Kjærgaard passed away in 2013. We also thank the leadership team at the Juliane Marie Center in Denmark for their support. The Rheumatology departments at the following hospitals in Denmark facilitated collection of data and samples: Rigshospitalet (Glostrup), Odense Universitetshospital, Dansk Gigthospital (Sønderborg), Aarhus University Hospital (Skejby), and Regionshospitalet Viborg. We thank all members of our project team for making this work possible: Anne-Grethe Rasmussen, Charlotte Schön Frengler, Dorte Heide, Randi Petersen, Tove Thorup Rasmussen, Lone Thomasen, Britta Hvidberg Nielsen, Teresa Rozenfeld, Kirsten Junker, Lis Kastberg Schubert, Lis Lund, Jette Barlach, Helle Bendtsen, Helle Andersen, and Marjo Westerdahl for their contribution with data and sample collection; Rikke Godtkjær Andersen, Mie Rams Rasmussen, Katrine Elmgaard Jensen, Pia Pedersen, Stine Birkelund, Louise Mielke, and Andreas Smed for management of data and samples. We also greatly appreciate valuable assistance provided by DANBIO personnel. DJ accepts responsibility for the integrity and validity of the data collected and analyzed.
Footnotes
This work was supported in part by funds from the National Institutes of Arthritis, Musculoskeletal and Skin Diseases, USA (Grants R21AR057931 and R01AR073111); Gigtforeningen, Denmark (Grant R87-A1477-B512); the Juliane Marie Center, Rigshospitalet, Denmark; and a private donor. These funders did not have any role in conducting this study or in the interpretation and reporting of results.
AP and MW contributed equally to this work.
The authors declare no conflicts of interest.
- Accepted for publication December 3, 2020.
- Copyright © 2021 by the Journal of Rheumatology