Peripheral Blood RNA Expression Profiling in Illicit Methcathinone Users Reveals Effect on Immune System

Methcathinone (ephedrone) is relatively easily accessible for abuse. Its users develop an extrapyramidal syndrome and it is not known if this is caused by methcathinone itself, by side-ingredients (manganese), or both. In the present study we aimed to clarify molecular mechanisms underlying this condition. We used microarrays to analyze whole-genome gene expression patterns of peripheral blood from 20 methcathinone users and 20 matched controls. Gene expression profile data were analyzed by Bayesian modeling and functional annotation. Of 28,869 genes on the microarrays, 326 showed statistically significant differential expression with FDR adjusted p-values below 0.05. Quantitative real-time PCR confirmed differential expression for the most of the genes selected for validation. Functional annotation and network analysis indicated activation of a gene network that included immunological disease, cellular movement, and cardiovascular disease functions (enrichment score 42). As HIV and HCV infections were confounding factors, we performed additional stratification of subjects. A similar functional activation of the “immunological disease” category was evident when we compared subjects according to injection status (past versus current users, balanced for HIV and HCV infection). However, this difference was not large therefore the major effect was related to the HIV status of the subjects. Mn–methcathinone abusers have blood RNA expression patterns that mostly reflect their HIV and HCV infections.


INTRODUCTION
Methcathinone (mcat) is an illicit stimulant drug synthesized by oxidation of ephedrine or pseudoephedrine with potassium permanganate. Motor dysfunction, such as hypokinesia and impaired balance, in methcathinone users has been reported in several case series from different countries, and probably develops due to chronic manganese intoxication (Sikk et al., 2007;Selikhova et al., 2008;Stepens et al., 2008). As methcathinone is relatively easy to manufacture using household materials, this poses a significant public health risk. Methcathinone abuse and its preparation with potassium permanganate is quite common in eastern European countries, but there are already signs that the practice is spreading to western countries (De Bie et al., 2007).
Molecular mechanisms underlying combined Mn and methcathinone toxicity are poorly understood. Mn neurotoxicity is well known and proposed mechanisms for neuronal death are described in several review articles (Fitsanakis et al., 2006;Burton and Guilarte, 2009;Curran et al., 2009). However, methcathinone has a molecular structure similar to methamphetamine and its pharmacodynamic properties suggest potential toxicity as well (Glennon et al., 1995;Rockhold et al., 1997;McCann et al., 1998). There are no molecular studies on the effect of these substances in combination.
Gene expression profiling using whole-genome microarray screening is a well-established technology to obtain a snapshot of the complex response variables (transcriptome) in biological systems. Expression profiling has great potential for examining molecular pathogenesis and biomarkers of disease. It is a relatively easy to use technology which gives a genome-wide overview of the molecular changes during pathological conditions. Peripheral blood gene expression has been used as a surrogate fingerprint of cerebral neurological diseases (Sharp et al., 2006;Sullivan et al., 2006;Mohr and Liew, 2007;Davies et al., 2009). The resulting profile is a complex description of the molecular phenotype and could be used to analyze different acute and chronic conditions.
The aim of this study was to analyze gene expression profiles in chronic methcathinone abusers. In addition, we applied functional genomic analysis to identify gene networks related to this chronic condition.

STUDY PARTICIPANTS
Twenty illicit methcathinone users and 20 sex-and age-matched controls reporting to be non-drug users and not to be infected with HIV and HCV were enrolled in the present study. Methcathinone was the main drug of abuse. Subjects were divided according to self reported drug abuse history into current or past users (reported having discontinued at least 1 year previously). Main subject characteristics including the results of motor and mental function and quality of life assessment scales are summarized in Table 1 (Hoehn and Yahr, 1967;Schwab and England, 1969;Folstein et al., 1975;Fahn et al., 1987;Peto et al., 1995). All subjects were positive for hepatitis C virus (HCV). Written informed consent was obtained from all subjects. The study was approved by the ethics committee on human research at Tartu University.

SAMPLE COLLECTION AND RNA PREPARATION
Blood samples from 20 subjects and 20 healthy controls were collected into Tempus tubes (Applied Biosystems, Foster City, USA). Blood was frozen and stored until further processing. RNA extraction from whole blood was performed according to the manufacturer's protocol (PN 4379228C). After RNA extraction, alpha and beta globin mRNA was depleted with GlobinClear Whole Blood Globin Reduction kit (Ambion, Austin, USA). The quality of RNA was analyzed with a Bioanalyzer 2100 (Agilent, Santa Clara, USA) and gene expression profiling was performed with GeneChip Human Gene 1.0 ST Arrays (Affymetrix, Santa Clara, USA) containing probes for all exons of 28,869 genes. We analyzed all the genes on the array without any filtering for presence/absence cells.

MICROARRAY HYBRIDIZATION AND ANALYSIS
In order to label the RNA we used the Affymetrix GeneChip Whole Transcript (WT) Sense Target Labeling Assay (Affymetrix, Santa Clara, USA), which is designed to generate amplified and biotinylated sense-strand targets from the entire expressed genome without bias. Briefly, double-stranded cDNA was synthesized from 300 ng of total RNA by reverse transcription using random hexamers tagged with a T7 promotor primer sequence. The doublestranded cDNA was subsequently used as a template and amplified by T7 RNA polymerase producing many copies of antisense cRNA. In the second cycle of cDNA synthesis, random hexamers were used to prime reverse transcription of the cRNA from the first cycle to produce single-stranded DNA in the sense orientation. This DNA was fragmented with a combination of uracil DNA glycosylase (UDG) and apurinic/apyrimidinic endonuclease 1 (APE 1). DNA was labeled by terminal deoxynucleotidyl transferase (TdT) and hybridization was performed according to the manufacturer's protocol. The arrays were subsequently washed, stained with phycoerythrin streptavidin and scanned according to standard Affymetrix protocol. Images were processed using the  Fahn et al., 1987). Hoehn and Yahr stages range from 1 to 5 (higher scores indicating more severe disability; Hoehn and Yahr, 1967). Schwab and England scores range from 0 to 100 (lower score indicating more severe disability; Schwab and England, 1969). PDQ-39 scores from 0 to 100 (higher score indicating worse quality of life; Peto et al., 1995). MMSE scores range from 0 to 30 (scores under 23 indicating dementia; Folstein et al., 1975).
*Denotes the subject without clinically relevant extrapyramidal syndrome.

Frontiers in Genetics | Neurogenomics
August 2011 | Volume 2 | Article 42 | 2 Affymetrix Expression Console and the MAS 5.0 algorithm was used to measure quality control parameters.

STATISTICAL ANALYSIS
The normalized, background subtracted, and modeled expression (Robust Multi-array Analysis, RMA) data (GEO accession number GSE28686) were further analyzed using a linear model combined with Bayesian moderation for standard errors implemented in the Bioconductor limma package of the statistical software R 1 (Ihaka and Gentleman, 1996;Smyth, 2004). False discovery rate 1 http://www.r-project.org/ (FDR) was used to address multiple testing corrections (Benjamini and Hochberg, 1995;Storey and Tibshirani, 2003). Moderated model calculates B-value, which is the log-odd that the gene is differentially expressed. B-statistics of zero corresponds to a 50-50 chance that the gene is differentially expressed and B-statistics is automatically adjusted for multiple testing.

QUANTITATIVE REAL-TIME PCR ANALYSIS
Eight genes with significant differences and with potentially interesting functions from the gene expression profiling data were further analyzed by means of real-time PCR (RT-PCR). These genes were: C15orf26, GPR15, IGFR1, SNRPN, IFI44L, IFI44, IFI27, and  (Smyth, 2004). For significance analysis moderated t-statistics combined with simple Bayesian model was used. "Probeset ID" is the Affymetrix ID for probes on the array, "Gene Symbol" is the official HUGO symbol for genes, "logFC" is log2 fold change, "AveExpr" is average expression value (log2), "adj. p-value" is FDR adjusted p-value, "B-value" is log-odds that the gene is differentially expressed.
www.frontiersin.org IFNG. RNA was converted into cDNA using the High Capacity cDNA Synthesis kit from Applied Biosystems. TaqMan assays and Gene Expression Master mix was used for the RT-PCR reaction in the SDS 7900 HT system (Applied Biosystems). Sample (20 controls and 20 methcathinone users in four replicates) comparisons were made using Welch's t -test.

FUNCTIONAL ANNOTATION OF DIFFERENTIALLY EXPRESSED GENES
To define the functional networks of the differentially expressed genes, data were analyzed by ingenuity pathway analysis (IPA, Ingenuity Systems 2 ). A data set containing Affymetrix probeset identifiers and corresponding fold change (log2) values was uploaded. Each gene identifier was mapped to its corresponding gene object in the Ingenuity Pathways Knowledge Base to generate 2 www.ingenuity.com the list of focus genes. Networks of these focus genes were then algorithmically generated based on inter-gene relationships manually curated in the Ingenuity database from published studies. Ingenuity pathway analysis calculates a significance score for each network. The score is generated using a p-value calculation, and is displayed as the negative logarithm of that p-value. This score indicates the likelihood that the assembly of a set of focus genes into a network could be explained by chance alone (e.g. score of 2 indicates that there is a 1˚in 100 chance that the focus genes are together in a network by random chance).

GENE EXPRESSION PROFILING WITH MICROARRAYS
Comparison of the blood RNA samples isolated from methcathinone users and healthy controls revealed distinct gene expression profiles. Statistically significant differential expression was observed for 326 genes with p-values less than 0.05 FIGURE 1 | Methcathinone users group distinctly on the hierarchical clustering heatmap of blood gene expression levels.
The top 100 genes from the decreasing ordered list of moderated t -values were clustered according to similarity in their gene expression patterns. Signals are scaled to Z -scores of the rows. The colored bar above the heatmap indicates the grouping variable -goldenrod for methcathinone user, blue for controls. Column labels at the bottom combine HIV status (positive or negative) with the drug user status (current, past, never). Row labels are gene symbols and the asterisk denotes the subject without clinically relevant extrapyramidal syndrome. Note some overlap between the two groups at the middle of the heatmap.
(FDR adjusted for multiple testing, Table 2). The magnitudes of expression signal differences (fold change or logFC) between these two groups were generally low -17 genes were upregulated in the methcathinone user more than 1.1-fold and only one gene was down-regulated more than 1.1-fold. However, the B-statistics values for the list of genes were very high suggesting real biological differences between these groups. A heatmap (Figure 1) from unsupervised hierarchical clustering indicates clear clustering of samples into drug users and non-users. However, in the middle of the heatmap, some overlap and mixture between groups occurs. These mcat samples are persons with HIV negative status. Therefore, HIV status has a significant impact on the gene expression profile and is a strong confounding factor. A functional annotation analysis was performed to define gene networks within the observed expression profiles. Relationships among the 326 genes significantly different in methcathinone users compared to control groups included two major networks. The first (with the highest score of 47; score indicates the −log of enrichment p-value) contained genes annotated with cell death, antigen presentation, or neurological disease functions (Table 3, Figure 2). A second network (score 42) was related to immunological disease, cellular movement, or cardiovascular disease. Enrichment or activation of these networks suggests that methcathinone users have changes in the genetic pathways related to the function of the nervous and immune systems. Table 2 includes many genes that are related to immune response, and may be differentially expressed compared to healthy controls because most of drug users were HIV positive. We therefore compared the subgroup of HIV negative drug users to the healthy control group. In this comparison, no genes were significantly different at FDR less than 5%, but 93 genes were different at or below a 30% FDR threshold ( Table A1 in Appendix). Annotation analysis for this dataset of HIV negative drug users indicated up-regulation (23 genes, score 54) in a network of genes with functions including hematological disease, cellular assembly, and organization or cell cycle ( Table 3). In another subgroup analysis we compared HIV positive subjects according to their user status, current versus discontinued. This comparison should indicate the difference caused by active drug abuse ( Table 3). In this case both groups are HIV and HCV positive, so infection status is balanced. The most significantly enriched network in this comparison included genes in the Genetic Disorder, Immunological Disease, or Cellular Movement categories (25 genes, enrichment score 46). This network supports the involvement of the immune system in the drug abuse induced pathologies.

QUANTITATIVE REAL-TIME PCR
To further verify changes found with the GeneChip experiment, we performed quantitative RT-PCR and analyzed gene expression levels of eight selected genes. For seven genes, statistically significant differential expression between 20 healthy controls and 20 methcathinone users was seen in the expected direction of change (C15orf26 p < 0.0001, GPR15 p < 0.0001, IGFR1 p = 0.0002, SNRPN p = 0.038, IFI44L p = 0.0264, IFI27  p = 0.0450, IFNG p < 0.0001). These results support the validity of our microarray gene expression analysis.

DISCUSSION
In this study, we examined changes in gene expression profiles in peripheral blood induced by methcathinone abuse. All but one of the subjects had a clinically relevant extrapyramidal syndrome. RNA expression patterns of drug users clearly differed from those of healthy controls. Also, RT-PCR results confirmed differences between these two groups. The overlap between controls and HIV negative methcathinone users in the middle of the heatmap (Figure 1) suggested that most of the hierarchical clustering was driven by HIV infection status. During recent years several microarray studies of HIV induced alterations in host cell gene expression have been performed using in vitro infected primary cells, cell lines, or tissue samples (Giri et al., 2006). However, few in vivo studies in HIV infected humans have evaluated expression profiles of tissue samples (Masliah et al., 2004;Everall et al., 2005) Frontiers in Genetics | Neurogenomics or peripheral blood mononuclear cells (Motomura et al., 2004;Kottilil et al., 2009;Monaco et al., 2009). Several of the same apoptosis, antigen presentation, and immune response genes (IFIT3,  IFI27, IFI44, IFNG, NFAT, STAT1, TCR, CCR5, ERK, IFITM3,  LY6E, IFNG, NFkB, OAS1, PLSCR1) described in previous HIV microarray studies were upregulated in our analysis. Thus peripheral blood can be used for finding gene expression alterations in HIV infected humans. Besides cell death and antigen presentation the genetic network with the highest score was also related to neurological disease, consistent with known associations of HIV infection and neurodegeneration (Clifford et al., 2005). Further minor subgroup analysis of HIV negative users who were all positive for HCV indicated that genes related to cell cycle, cellular growth and proliferation, cancer, and cellular development were upregulated compared to healthy controls. The same pathways have been shown to be involved in HCV induced hepatic changes like cirrhosis, dysplasia, or hepatocellular carcinoma (Wurmbach et al., 2007;De Giorgi et al., 2009).
As the aim of our study was to analyze the drug abuse effect, we decided to stratify our group of subjects. We analyzed only HIV and HCV positive subjects according to their injection status (past versus current) and a clinical measure of the severity of the neurological syndrome UPDRS. The most significantly enriched network after comparing HIV positive current and past users included genes annotated with functions in Genetic Disorder, Immunological Disease, or Cellular Movement categories (25 genes, enrichment score 46). This network confirms the involvement of the immune system in drug abuse induced pathologies. However, as both comparison groups had neurological symptoms we can not draw any conclusions on the potential causes of the extrapyramidal syndrome. On the other hand, injection status seems to influence immune system and this is one factor involved in the development of a neurodegenerative syndrome. Psychostimulant abuse has been shown to compromise immunological status (Everall et al., 2005). Chronic Mn exposure has similar effects on the immune system (Sengupta et al., 2007). Thus, both of these factors have potential to induce neurodegeneration.
The UPDRS score reflects the clinical severity of the extrapyramidal syndrome. No significant correlations were found between the duration of Mn-methcathinone use and UPDRS or between blood RNA profiles and UPDRS.
There are several limitations in our study. The main problem is related to the sample organization. Our sample contains many confounding factors -HIV and HCV infection, different duration of methcathinone use, concomitant consumption of other drugs, alcohol, and tobacco. We have not analyzed the clinical status of the HIV and HCV infections. In addition, several drug abusers were undernourished which makes identification of a relevant control group even more complicated. To overcome these limitations, we have applied complex bioinformatics analytical tools in relatively large samples (N = 20 in both groups) that allowed us to perform some cohort stratification. For further studies within methcathinone users more balanced control samples are needed. We have tried to collect similar balanced sample as a control group, but due to the specific and complex nature of our study group it is difficult and time consuming. Ideal controls would be either drug-free HIV and HCV positive subjects or addicts abusing other drugs. We fully understand the limitations caused by the complexity of the problem we analyze. On the other hand, the situation of a non-balanced control sample is quite common in clinical research. "Probeset ID" is the Affymetrix ID for probes on the array, "Gene Symbol" is the official HUGO symbol for genes, "logFC" is log2 fold change, "AveExpr" is average expression value (log2), "adj. p-value" is FDR adjusted p-value, "B-value" is log-odds that the gene is differentially expressed.