Article

Original Article

Exp Neurobiol 2016; 25(1): 40-47

Published online February 29, 2016

https://doi.org/10.5607/en.2016.25.1.40

© The Korean Society for Brain and Neural Sciences

Blood Transcriptome Profiling in Myasthenia Gravis Patients to Assess Disease Activity: A Pilot RNA-seq Study

Kee Hong Park1, Junghee Jung2, Jung-Hee Lee3 and Yoon-Ho Hong4*

1Department of Neurology, Gyeongsang National University Hospital, Jinju 52727, 2Department of Bioinformatics, Macrogen Inc., Seoul 08511, 3Department of Biomedical Science, Hallym University, Chuncheon 24252, 4Department of Neurology, Seoul Metropolitan Government Seoul National University Boramae Medical Center, Seoul 07061, Korea

Correspondence to: *To whom correspondence should be addressed.
TEL: 82-2-840-2474, FAX: 82-2-831-2826
e-mail: nrhong@gmail.com

Received: January 6, 2016; Revised: February 5, 2016; Accepted: February 5, 2016

Myasthenia gravis (MG) is an antibody-mediated autoimmune disease characterized by exertional weakness. There is no biomarker to reflect disease activity and guide treatment decision. Here, we reported a pilot blood transcriptome study using RNA sequencing (RNA-seq) that identified differences of 5 samples in active status and 5 in remission from 8 different patients and 2 patients provided samples for both active and remission phase. We found a total of 28 differentially expressed genes (DEGs) possibly related to disease activity (23 up-regulated and 5 down-regulated). The DEGs were enriched for the cell motion and cell migration processes in which included were ICAM1, CCL3, S100P and GAB2. The apoptosis and cell death pathway was also significantly enriched, which includes NFKBIA, ZC3H12A, TNFAIP3, and PPP1R15A. Our result suggests that transcript abundance profiles of the genes involved in cell trafficking and apoptosis may be a molecular signature of the disease activity in MG patients.

Keywords: Transcriptome, RNA sequencing, myasthenia gravis, cell migration, apoptosis

Myasthenia gravis (MG) is an antibody-mediated autoimmune disease characterized by exertional skeletal muscle weakness, caused by auto-antibodies to components of muscle membrane at neuromuscular junction (NMJ) [1]. Acetylcholine receptor (AChR) antibodies are detected in about 85% of MG patients, and antibodies to muscle-specific tyrosine kinase (MuSK) are positive in approximately half of the generalized seronegative MG patients [2]. The list of pathogenic auto-antibodies in MG has been recently expanded with discovery of antibodies to low-density lipoprotein receptor-related protein 4 (LRP4) and agrin in some of the double-seronegative MG patients [3]. These highly specific auto-antibodies in MG have a crucial role for the diagnosis of disease, and MuSK antibody in particular seem to be useful in predicting adverse effects to acetylcholinesterase inhibitors and refractoriness to conventional immunosuppressive treatment [4]. Although the levels of these auto-antibodies tend to decrease with immunosuppressive treatment, they are highly variable between patients, and do not correlate with disease severity [5]. There is an unmet clinical need to develop biomarkers that can reflect disease activity and/or severity, and to guide treatment decisions.

Gene expression profiling on a genome-wide scale has been increasingly used to investigate pathogenesis, and also to develop potential biomarkers in various human diseases. Advances of high-throughput techniques such as DNA microarrays and recently RNA sequencing (RNA-seq) allow us to analyze gene expression in an unbiased and comprehensive way [6]. In MG, distinct gene expression signatures were found in the thymus of AChR antibody positive patients using microarrays, which include dysregulation of chemokines such as CCL21, CXCL13, CXCL10 and CXCR3 [7,8,9]. Interferon type I overexpression together with dysregulated expression of dsRNA signaling molecules was found in thymoma-associated MG patients [10]. Of note, previous gene expression studies mostly have used thymic tissues which are not readily accessible. Furthermore, it is not feasible to sample repeatedly in order to monitor disease activity and response to treatment. Since the status of human immune system may be best monitored by the changes of the composition and transcript abundance of circulating immune cells, blood transcriptome could be a good alternative target of investigation. Indeed, recent studies in autoimmune disease such as systemic lupus erythematosus, multiple sclerosis, and psoriasis have demonstrated perturbation of blood transcriptome, and identified molecular signatures to predict clinical relapse and response to treatment [11,12,13]. Here, we report the results of blood whole transcriptome study using RNA-seq in 8 patients with MG who were in either remission or active disease states. We found 28 genes that are differentially expressed (DEGs) according to disease activity. Functional analysis suggested that genes expression profiles related to immune cell trafficking and apoptosis might be a molecular signature of the disease activity in MG.

Patients

This study was approved by the Institutional Review Board of the Seoul Metropolitan Government Boramae Medical Center (IRB no. 16-2013-116). Written informed consent was obtained from the participants who were diagnosed as having MG and on the status of either active or remission. The active group consisted of 5 patients with de novo or refractory MG (moderate to severe symptoms despite long-term immunosuppressive treatment). Disease severity was graded according to the Myasthenia Gravis Foundation of America (MGFA) Clinical Classification [14]. Remission was defined by the MGFA post-intervention status, and included complete stable remission (CSR), pharmacologic remission (PR), and minimal manifestation (MM) (Table 1). Two patients provided samples at different time points, one during active disease status and the other in remission state. There were not statistically differences of mean age (p=0.69), disease duration (p=0.31), and AChR antibody titer (p=0.69) between 2 groups.

PBMC isolation and RNA purification

For isolation of peripheral blood mononuclear cells (PBMC), the Lymphoprep™ was used according to the manufacturer's protocol (Axis-shield, Oslo, Norway). Medium was placed in the tube, and then blood sample diluted with saline with 1:1 was added. After centrifugation for 20 minutes at 600×g , sedimented PBMCs were harvested. RNA purification was performed with the RNeasy Mini kit with the isolated PBMC sample (Qiagen, Seoul, Korea). The cell pellet was mixed with RLT buffer and 70% ethanol. The lysate was then loaded onto the RNeasy Mini spin column to facilitate the binding of RNA to the column and for the removal of contaminants. DNase was added to remove residual DNA efficiently.

RNA-Seq

The mRNA-Seq sample was obtained using Illumina TruSeq™ RNA Sample Preparation Kit (Illumina, Inc., San Diego, CA, USA). In brief, purifying the poly-A containing mRNA molecules with poly-T oligo-attached magnetic beads was the first step, followed by thermal mRNA fragmentation. The RNA fragments were then transcribed into first strand cDNA using reverse transcriptase and random primers. The cDNA was synthesized to second strand cDNA using DNA Polymerase I and RNase H. After the end repair process, single 'A' bases were added to the fragments and adapters were then ligated, preparing cDNA for hybridization onto a flow cell. Finally, the products were purified and enriched with PCR to create the cDNA library (Macrogen, Seoul, Korea).

Aligning RNA-Seq reads and abundance estimation

Fragmented cDNAs were aligned using TopHat v.2.0.11 [15] and subsequently aligned with sequences obtained from the human genome (UCSC version hg19) using the Bowtie 2.1.0 algorithm [16]. Abundance of aligned reads were estimated by Cufflinks v.2.1.1 [17], which accepted aligned reads and assembled the alignments into a simple and clear set of transcripts. Next, RNA-seq fragment counts were measured by the unit of fragments per kilobase of exon per million fragments mapped (FPKM) [18]. DESeq, another tool for DEG analysis, was used to compare the results with Cuffdiff analysis. Cuffdiff determines differential expression using t-test from FPKM values and is based on beta negative binomial model [19], while DESeq uses exact test based on negative binomial model [20]. We compared the results from Cuffdiff and DESeq analyses, and took the intersection of them for downstream pathway analysis.

Statistical analysis

For DEG analysis, the values of log2 (FPKM+1) were calculated, and these were normalized by quantile normalization. p-values were obtained by t-test between the active and remission groups, and fold changes were calculated with the mean log2 (FPKM+1) values, gene by gene. All data analysis of DEG was conducted using R 2.14.1 (http://www.r-proj ect.org). To segregate the samples according to the disease activity, a multi-dimensional scaling (MDS) analysis was done.

Pathway analysis using DAVID and IPA

For functional enrichment analysis using gene ontology (GO), the Database for Annotation, Visualization and Integrated Discovery (DAVID v.6.7) was used. The list of commonly detected genes both in Cuffdiff and DESeq analysis was uploaded via the web interface (http://david.abcc.ncifcrf.gov), and the background was designated as Homo sapiens [21]. Functional annotation clusters were selected using the cutoff threshold enrichment score 1.3 which is equivalent to a non-log scale score of 0.05. Gene interaction was visualized by Gene Multiple Association Network Integration Algorithm (GeneMANIA) (http://www.genemania.org/). Another functional analysis was performed with the Ingenuity Pathway Analysis (IPA) software (licensed use of Ingenuity Systems, http://www.qiagen.com/ingenuity).

Differentially expressed genes (DEG)

In total, there were 48,385 transcripts, and we excluded any transcripts with an FPKM value of '0,' leaving 10,640 transcripts to be analyzed. MDS analysis showed that our MG samples could be differentiated according to the disease activity along the first plot dimension (Fig. 1).

The level of expression of 98 genes was significantly different between the two groups (fold change≥2, p-value<0.05) with 63 genes up-regulated and 35 genes down-regulated in the remission group (Supplementary Table 1, Fig. 2). The number of DEGs derived with DESeq analysis was greater with 165 up-regulated and 127 down-regulated in the remission group (Supplementary Table 2). Twenty-eight genes were common in the results of both Cuffdiff and DEseq analyses (23 up-regulated and 5 down-regulated genes, Table 2, Fig. 3). When only the two paired samples were analyzed separately, we found 18 up-regulated genes (fold change≥2, p-value<0.05), with no down-regulated genes (Supplementary Table 3). All the up-regulated genes from the paired sample analysis except for TPGS1 were contained in the list of DEGs analyzed with the whole study samples.

Pathway analysis using DAVID and IPA

DAVID analysis with up-regulated genes revealed two enriched GO functional annotation clusters (Table 3). The first cluster consisted of biological functions of cell motion and cell migration. The second cluster consisted of apoptosis and cell death. Number of down-regulated genes was not sufficient to analyze. Gene interaction of up-regulated genes was shown in Fig. 4.

According to the downstream effect analysis performed by using IPA, top 5 functions affected were hematologic system development and function, immune cell trafficking, cellular movement, inflammatory response and cell-to-cell signaling and interaction pathways. Most of the genes related to these categories were down-regulated in active group.

In the present study, whole transcriptome RNA sequencing was performed in 10 PBMC samples obtained from eight MG patients to investigate systemic changes of gene expression possibly relevant to disease activity. Because there is no general consensus regarding the best method for the differential expression analysis of RNA-seq data, two different analysis tools (Cuffdiff and DESeq) were used separately with the intersection of the results being taken for downstream pathway analysis. The results demonstrated that among 10,640 transcripts only 28 genes were differentially expressed with 23 being up-regulated and 5 down-regulated in remission compared to active status of the disease. Interestingly, major functional themes of the differentially expressed genes include differentiation, trafficking and apoptosis of immune cells, inflammatory response and cell-to-cell signaling.

S100B (S100 calcium-binding protein B), which was the most down-regulated gene in the remission group, is a ligand of the receptor for advanced glycation end products (RAGE). In the experimental autoimmune myasthenia gravis (EAMG) model, S100B levels were significantly higher [22]. Mean value of S100B of MG was higher than the normal controls, but it was not statistically significant [23]. CTTN (cortactin) was another down-regulated gene in the remission group, which is needed for the formation of the AChR cluster and antibody against cortactin was reported in MG patients [24]. ABL1 (Abelson murine leukemia viral oncogene homolog 1) is up-regulated in the remission group, which is critical mediator of postsynaptic assembly at the NMJ via providing a specific tyrosine kinase activity downstream of the MuSK receptor that is required for agrin-induced AChR clustering [25].

One of the significant functional annotation cluster in DAVID analysis was the apoptosis and programmed cell death. Immunosuppressive agents are the mainstay of treatment for MG and their mechanisms are related with these pathways. Corticosteroids induce T lymphocyte apoptosis, and other agents such as azathioprine also inhibit cell proliferation [26]. All of the remission group patients received prednisolone and three patients also received tacrolimus, another T cell apoptosis inducer [27]. CXCL13 and CCL21 were overexpressed in hyperplastic thymus [8,28] and CXCR5 was overexpressed on T cells of MG patients [29]. Most of the genes associated with chemotaxis, however, were down-regulated in the active group of the present study. This discrepancy needs to be confirmed in further studies, but might be explained by the followings. First, the previous studies had compared between MG patients and healthy controls. Second, difference of the target samples (thymus or neuromuscular junction in the previous studies) might account for the discrepancy. In a previous study of rheumatoid arthritis, chemokine receptor profile of PBMCs was increased with treatment, which the authors explained as systemic compensation for the changes in the inflamed tissue [30].

There are several limitations to be acknowledged in this pilot study. Sample size was too small to make a robust conclusion. Immunosuppressive treatment varied among patients, which might affect the results. Despite these limitations, to our knowledge, this is the first study investigating the systemic changes of blood whole transcriptome according to the disease activity in MG patients. Genes that are involved in the biological process of immune cell trafficking and apoptosis might be promising candidates for biomarkers, and warrant further studies.

Supplementary Table 1

List of differentially expressed genes derived with Cuffdiff.

en-25-40-s001.pdf

Supplementary Table 2

List of differentially expressed genes derived with DESeq.

en-25-40-s002.pdf

Supplementary Table 3

List of up-regulated genes in the remission group of the paired samples.

en-25-40-s003.pdf
Fig. 1. Multidimensional Scaling plot of gene expression profiles. The first plot dimension roughly corresponds to the disease activity. Paired samples are linked by dotted arrows.
Fig. 2. Differentially expressed genes by Cuffdiff. (A) Volcano plot illustrating the differential expression levels of genes of the active and remission group. Genes that are significantly up- and down-regulated in the active compared to remission group are shown in red and green dots, respectively. (B) Heat map of the hierarchical clustering based on 98 differentially expressed genes (fold change≥2, p-value<0.05).
Fig. 3. Venn diagram of DEGs from DESeq and Cuffdiff. (A) Up-regulated genes in the remission group. (B) Down-regulated genes in the remission group.
Fig. 4. Gene interactions of up-regulated genes from GeneMANIA. Blue line indicates co-localization of the genes and predicted functional relationships between genes are indicated in orange. Circles filled with black signify the commonly detected genes by Cuffdiff and DESeq. Circles filled with gray indicate their interaction and added by GeneMania.
Table. 1. Demographics of study population

*Patient 1 and 6 are same patients..

Patient 2 and 7 are same patients..

Abbreviations: MG: myasthenia gravis; AChR Ab: acetylcholine receptor antibody; MGFA: Myasthenia Gravis Foundation of America; M: male; F: female; NL: normal; Pd: prednisolone; MM: minimal manifestation; PR: pharmacologic remission..


Table. 2. List of differentially expressed genes (from both Cuffdiff and DESeq analyses)
Table. 3. Functional annotation cluster analysis using DAVID

*Gene Ontology..


  1. Berrih-Aknin S, Le Panse R. Myasthenia gravis: a comprehensive review of immune dysregulation and etiological mechanisms. J Autoimmun 2014;52:90-100.
    Pubmed
  2. McConville J, Farrugia ME, Beeson D, Kishore U, Metcalfe R, Newsom-Davis J, Vincent A. Detection and characterization of MuSK antibodies in seronegative myasthenia gravis. Ann Neurol 2004;55:580-584.
    Pubmed
  3. Berrih-Aknin S, Frenkian-Cuvelier M, Eymard B. Diagnostic and clinical classification of autoimmune myasthenia gravis. J Autoimmun 2014;48-49:143-148.
    Pubmed
  4. Evoli A, Tonali PA, Padua L, Monaco ML, Scuderi F, Batocchi AP, Marino M, Bartoccioni E. Clinical correlates with anti-MuSK antibodies in generalized seronegative myasthenia gravis. Brain 2003;126:2304-2311.
    Pubmed
  5. Heldal AT, Eide GE, Romi F, Owe JF, Gilhus NE. Repeated acetylcholine receptor antibody-concentrations and association to clinical myasthenia gravis development. PLoS One 2014;9:e114060.
    Pubmed
  6. Chaussabel D, Pascual V, Banchereau J. Assessing the human immune system through blood transcriptomics. BMC Biol 2010;8:84.
    Pubmed
  7. Le Panse R, Cizeron-Clairac G, Bismuth J, Berrih-Aknin S. Microarrays reveal distinct gene signatures in the thymus of seropositive and seronegative myasthenia gravis patients and the role of CC chemokine ligand 21 in thymic hyperplasia. J Immunol 2006;177:7868-7879.
    Pubmed
  8. Meraouna A, Cizeron-Clairac G, Panse RL, Bismuth J, Truffault F, Tallaksen C, Berrih-Aknin S. The chemokine CXCL13 is a key molecule in autoimmune myasthenia gravis. Blood 2006;108:432-440.
    Pubmed
  9. Feferman T, Aricha R, Menon R, Souroujon MC, Berrih-Aknin S, Fuchs S. DNA microarray in search of new drug targets for myasthenia gravis. Ann N Y Acad Sci 2007;1107:111-117.
    Pubmed
  10. Cufi P, Soussan P, Truffault F, Fetouchi R, Robinet M, Fadel E, Berrih-Aknin S, Le Panse R. Thymoma-associated myasthenia gravis: on the search for a pathogen signature. J Autoimmun 2014;52:29-35.
    Pubmed
  11. Bauer JW, Petri M, Batliwalla FM, Koeuth T, Wilson J, Slattery C, Panoskaltsis-Mortari A, Gregersen PK, Behrens TW, Baechler EC. Interferon-regulated chemokines as biomarkers of systemic lupus erythematosus disease activity: a validation study. Arthritis Rheum 2009;60:3098-3107.
    Pubmed
  12. Hecker M, Goertsches RH, Fatum C, Koczan D, Thiesen HJ, Guthke R, Zettl UK. Network analysis of transcriptional regulation in response to intramuscular interferon-β-1a multiple sclerosis treatment. Pharmacogenomics J 2012;12:134-146.
    Pubmed
  13. Jabbari A, Suárez-Fariñas M, Fuentes-Duculan J, Gonzalez J, Cueto I, Franks AG, Krueger JG. Dominant Th1 and minimal Th17 skewing in discoid lupus revealed by transcriptomic comparison with psoriasis. J Invest Dermatol 2014;134:87-95.
    Pubmed
  14. Jaretzki A, Barohn RJ, Ernstoff RM, Kaminski HJ, Keesey JC, Penn AS, Sanders DB, Task Force of the Medical Scientific Advisory Board of the Myasthenia Gravis Foundation of America. Myasthenia gravis: recommendations for clinical research standards. Neurology 2000;55:16-23.
    Pubmed
  15. Trapnell C, Pachter L, Salzberg SL. TopHat: discovering splice junctions with RNA-Seq. Bioinformatics 2009;25:1105-1111.
    Pubmed
  16. Langmead B, Trapnell C, Pop M, Salzberg SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol 2009;10:R25.
    Pubmed
  17. Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, Salzberg SL, Wold BJ, Pachter L. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol 2010;28:511-515.
    Pubmed
  18. Twine NA, Janitz K, Wilkins MR, Janitz M. Whole transcriptome sequencing reveals gene expression and splicing differences in brain regions affected by Alzheimer's disease. PLoS One 2011;6:e16266.
    Pubmed
  19. Seyednasrollah F, Laiho A, Elo LL. Comparison of software packages for detecting differential expression in RNA-seq studies. Brief Bioinform 2015;16:59-70.
    Pubmed
  20. Anders S, Huber W. Differential expression analysis for sequence count data. Genome Biol 2010;11:R106.
    Pubmed
  21. Huang da W, Sherman BT, Lempicki RA. Systematic and integrative analysis of large gene lists using DAVID bioinformatics resources. Nat Protoc 2009;4:44-57.
    Pubmed
  22. Mu L, Zhang Y, Sun B, Wang J, Xie X, Li N, Zhang J, Kong Q, Liu Y, Han Z, Wang G, Fu Z, Yu B, Li G, Li H. Activation of the receptor for advanced glycation end products (RAGE) exacerbates experimental autoimmune myasthenia gravis symptoms. Clin Immunol 2011;141:36-48.
    Pubmed
  23. Moser B, Bekos C, Zimprich F, Nickl S, Klepetko W, Ankersmit J. The receptor for advanced glycation endproducts and its ligands in patients with myasthenia gravis. Biochem Biophys Res Commun 2012;420:96-101.
    Pubmed
  24. Gallardo E, Martínez-Hernández E, Titulaer MJ, Huijbers MG, Martínez MA, Ramos A, Querol L, Díaz-Manera J, Rojas-García R, Hayworth CR, Verschuuren JJ, Balice-Gordon R, Dalmau J, Illa I. Cortactin autoantibodies in myasthenia gravis. Autoimmun Rev 2014;13:1003-1007.
    Pubmed
  25. Finn AJ, Feng G, Pendergast AM. Postsynaptic requirement for Abl kinases in assembly of the neuromuscular junction. Nat Neurosci 2003;6:717-723.
    Pubmed
  26. Sanders DB, Evoli A. Immunosuppressive therapies in myasthenia gravis. Autoimmunity 2010;43:428-435.
    Pubmed
  27. Ponseti JM, Gamez J, Azem J, López-Cano M, Vilallonga R, Armengol M. Tacrolimus for myasthenia gravis: a clinical study of 212 patients. Ann N Y Acad Sci 2008;1132:254-263.
    Pubmed
  28. Berrih-Aknin S, Ruhlmann N, Bismuth J, Cizeron-Clairac G, Zelman E, Shachar I, Dartevelle P, de Rosbo NK, Le Panse R. CCL21 overexpressed on lymphatic vessels drives thymic hyperplasia in myasthenia. Ann Neurol 2009;66:521-531.
    Pubmed
  29. Saito R, Onodera H, Tago H, Suzuki Y, Shimizu M, Matsumura Y, Kondo T, Itoyama Y. Altered expression of chemokine receptor CXCR5 on T cells of myasthenia gravis patients. J Neuroimmunol 2005;170:172-178.
    Pubmed
  30. Nissinen R, Leirisalo-Repo M, Peltomaa R, Palosuo T, Vaarala O. Cytokine and chemokine receptor profile of peripheral blood mononuclear cells during treatment with infliximab in patients with active rheumatoid arthritis. Ann Rheum Dis 2004;63:681-687.
    Pubmed