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.

Materials and Methods

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.

Table 1

Clinical and demographic characteristics.

Subject no.	Sex	Age (years)	HIV	History of drug abuse	Duration of use (years)	UPDRS	Hoehn and Yahr	Schwab and England	PDQ-39	MMSE
1	M	47	pos	Current	6	13	2,5	95	24	27
2	F	24	pos	Past	1	27	2,5	80	23	30
3	M	44	pos	Current	10	36	3	80	57	28
4	M	43	pos	Current	4	17	2	90	9	25
5*	M	42	pos	Current	1	1	0	100	0	30
6	M	28	pos	Current	2	19	2	95	53	29
7	M	37	pos	Current	8	46	3	70	64	28
8	F	29	pos	Current	11	12	2	90	46	29
9	M	35	neg	Current	8	44	2	90	66	29
10	M	32	pos	Past	2	73	4	55	61	30
11	M	33	neg	Current	9	79	4	55	71	28
12	M	32	pos	Current	13	50	3	70	55	28
13	M	36	neg	Past	8	59	4	50	43	30
14	F	18	pos	Past	2	22	3	75	44	30
15	M	28	pos	Current	2	29	2,5	90	16	30
16	M	29	neg	Current	2	42	3	80	39	30
17	M	44	neg	Past	10	69	4	60	42	30
18	M	32	pos	Past	1	26	3	85	60	30
19	M	28	pos	Past	6	53	3	75	62	28
20	M	32	pos	Current	9	61	2,5	85	51	30

UPDRS scores range from 0 to 176 (higher scores indicating more severe disability; Fahn et al., .

*Denotes the subject without clinically relevant extrapyramidal syndrome.

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 double-stranded 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 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 (Ihaka and Gentleman, 1996; Smyth, 2004). False discovery rate (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 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). 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 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).

Results

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 (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.

Table 2

Differential expression analysis of the blood RNA suggests significant differences between study groups.

Probeset ID	Gene symbol	logFC	AveExpr	Adj. p-value	B-value	Gene title
8081214	GPR15	1.6	8.2	1.11E-04	10.199	G protein-coupled receptor 15
8046488	CDCA7	0.6	6.5	1.11E-04	9.901	Cell division cycle associated 7
7917283	MCOLN2	1.0	7.5	4.31E-04	8.329	Mucolipin 2
7964787	IFNG	1.1	7.1	0.001	7.697	Interferon, gamma
8081799	TIGIT	1.1	8.4	0.001	7.607	T cell immunoreceptor with Ig and ITIM domains
8038861	SIGLEC6	0.3	7.3	0.001	6.485	Sialic acid binding Ig-like lectin 6
8165653	ND1	0.4	13.2	0.001	6.463	NADH dehydrogenase, subunit 1 (complex I)
8175365	–	0.5	6.5	0.001	6.337	–
8165674	ND3	0.3	13.0	0.002	6.065	NADH dehydrogenase, subunit 3 (complex I)
8165648	C7orf11	0.5	12.5	0.002	6.027	Chromosome 7 open reading frame 11
7946563	–	0.6	11.3	0.002	5.917	–
7898793	C1QA	0.5	6.0	0.002	5.834	Complement component 1, q subcomponent, A chain
7903358	VCAM1	0.8	6.0	0.002	5.722	Vascular cell adhesion molecule 1
8149979	C8orf80	0.5	7.3	0.002	5.684	Chromosome 8 open reading frame 80
8169249	MID2	0.5	6.7	0.002	5.367	Midline 2
7983365	TRIM69	0.5	5.9	0.003	5.237	Tripartite motif-containing 69
8044154	CD8B	0.9	9.1	0.003	5.172	CD8b molecule
7944739	CRTAM	1.0	7.6	0.003	4.949	Cytotoxic and regulatory T cell molecule
8053584	CD8A	0.9	9.9	0.003	4.873	CD8a molecule
8109639	PTTG1	0.5	7.8	0.003	4.751	Pituitary tumor-transforming 1
7902660	WDR63	0.6	5.3	0.003	4.749	WD repeat domain 63
8019842	TYMS	0.8	6.7	0.003	4.682	Thymidylate synthetase
7948420	FABP5	1.0	6.5	0.003	4.662	Fatty acid binding protein 5 (psoriasis-associated)
8102643	CCNA2	0.6	6.8	0.003	4.653	Cyclin A2
7917576	GBP5	1.0	11.1	0.004	4.466	Guanylate binding protein 5
8002975	CDYL2	0.3	7.6	0.004	4.465	Chromodomain protein, Y-like 2
7930577	CASP7	0.5	8.1	0.004	4.455	Caspase 7, apoptosis-related cysteine peptidase
8151334	MSC	0.6	7.0	0.004	4.428	Musculin (activated B-cell factor-1)
7898535	–	0.3	3.9	0.004	4.367	–
8018352	SLC25A19	0.3	6.7	0.004	4.307	Solute carrier family 25 member 19
8053690	IGKC	1.0	9.4	0.004	4.291	Immunoglobulin kappa constant
7981962	SNRPN	0.6	7.1	0.004	4.159	Small nuclear ribonucleoprotein polypeptide N
7898805	C1QB	0.5	6.8	0.005	4.040	Complement component 1, q subcomponent, B chain
7940028	SERPING1	1.4	8.7	0.005	3.978	Serpin peptidase inhibitor, clade G (C1 inhibitor), member 1
7986359	IGF1R	0.7	8.6	0.005	3.914	Insulin-like growth factor-1 receptor
8116130	FAM153B	−0.6	6.9	0.005	3.898	Family with sequence similarity 153, member B
8043480	IGKV1OR15-118	0.9	10.9	0.005	3.863	Immunoglobulin kappa variable 1/OR15-118 pseudogene

Statistical analysis of Affymetrix GeneChip data was performed with the .

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.

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 3

Genetic networks significantly changed.

Molecules in network	Score	Focus molecules	Top functions
IN METHCATHINONE USERS
BAK1,CD3,CHI3L2,ENG, EOMES, FBXO6, GBP5, GBP6, GBP7, GBP4 (includes EG:115361), GCH1, growth hormone, HLA-DR, IFI27, IFI30, IFI44L, IFIT3, IFNG, KLF10, MT1F	47	28	Cell death, antigen presentation, neurological disease
AGER, C2, C1q, C1QA, C1QB, CCR5, CD80, CD81, CD40LG, complement component 1, Creb, CXCR6 (includes EG:10663), ERK,IFI44, IFITM3, IFNα/&beta	42	24	Immunological disease, cellular movement, cardiovascular disease
In HIV negative methcathinone users, compared to controls
Actin, adaptor protein 2, ANK1, Ap1, Ap2 alpha, AP2A1,calmodulin, caspase, CKS2, clathrin, CLU, DAB2, EGF, EIF2AK1, EPB41, F actin, FKBP8, GNAS, GUK1, HDGF, HSPB1, IgG, MPRIP, NFE2, NFkB (complex), NRGN	54	23	Hematological disease, cellular assembly, and organization or cell cycle
AGTR1B, ANKH, BAT3, betaestradiol, C14ORF45, CA2, CDKN1B, CKS2, CTSA, CUL2, CUL4A, DGKD (includes EG:8527), GRINA, GYPA, HNF4A, MIR135A1, MIRLET7C, MMD, NBL1, NEU1, PHOSPHO1, progesterone	31	15	Amino acid metabolism, cancer, cell morphology
In HIV positive methcathinone users, compared by user status
AICDA, CARD8, CCNH, CD40, CD226, CXCL5, DDR1, EIF2AK1, EMP1, ERN1, GYPA, IFIT1L, IFN beta, Iga, Ige, IgG, IGHA1, IGJ, IGLL1, Igm, IL1, IL12 (complex), IL1RAP, interferon alpha, IRAK1BP1, KLF3, MAP4K5, NFkB (complex), NLRP1, RIOK3, Sapk, THBD, TNFSF4, XCL1, ZNF675	46	25	Genetic disorder, immunological disease, cellular movement
14-3-3, adaptor protein 2, Ap1, ARHGEF12,CAMK2D, CCR3, CDC25A, CREG1, CSPG4, DAB2, DENND4A, DLEU2, E2f, ERK, GAB1, hCG, HDLBP, histone h3, ITGB3, LDLRAP1, LRRFIP1, Mapk, MAX, MYO10, ORC1L, PDGF BB, PGRMC1, Pkc(s), Ras, RGS13, RICS, SKP1, STON2, TFDP2, YBX1	43	24	Cellular assembly and organization, cancer, cell cycle

Figure 2

Annotation enrichment analysis indicated that a network including “cell death, antigen presentation, neurological disease” gene functions is significantly enriched in the subjects with methcathinone abuse [score 47, score is −log (. Red nodes designate upregulated genes and the number indicates log2 fold change (0 is for equal expression). Uncolored nodes are genes in this network that were not in our list of differentially expressed genes.

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.

Table A1

Expression analysis of the blood RNA between HIV negative drug users and the healthy control group.

Probeset ID	Gene symbol	logFC	AveExpr	adj. p-value	B-value	Gene title
7985352	C15orf26	−0.5	6.7	0.3736	−1.890	Chromosome 15 open reading frame 26
7982090	SNRPN	0.7	2.9	0.9245	−2.305	Small nucleolar RNA, C/D box 115-8
8018908	DNAH17	0.4	4.9	0.9245	−2.448	Dynein, axonemal, heavy chain 17
7942983	RNU6-2	0.3	4.0	0.9245	−2.574	RNA, U6 small nuclear 2
7893226	–	−0.9	6.1	0.9245	−2.588
7982068	SNRPN	0.6	2.9	0.9245	−2.784	Small nucleolar RNA, C/D box 115-8
8059578	–	0.6	6.6	0.9245	−2.791	ENSG00000208057
7982056	SNRPN	0.7	3.5	0.9245	−2.821	Small nucleolar RNA, C/D box 115-8
8123336	SNORD45A	0.4	4.5	0.9245	−2.888	Small nucleolar RNA, C/D box 45A
8167947	–	−1.3	4.4	0.9245	−2.894	ENSG00000211140
8081214	GPR15	1.1	7.6	0.9245	−2.922	G protein-coupled receptor 15
8146635	–	0.7	2.5	0.9245	−2.953	ENSG00000211143
7946054	HBG2	−1.9	11.1	0.9245	−2.963	Hemoglobin, gamma G
8165653	NDUFS1	−0.3	13.3	0.9245	−2.971	NADH-ubiquinone oxidoreductase chain 1
7922200	SELP	−0.7	8.5	0.9245	−2.975	Selectin P (granule membrane protein 140 kDa, antigen CD62)
8051571	–	0.4	4.4	0.9245	−3.019	ENSG00000211163
7961102	CLEC1B	−0.7	8.0	0.9245	−3.036	C-type lectin domain family 1, member B
7912582	LOC441873	0.6	5.1	0.9245	−3.048	PRAME family member 11
8015187	KRTAP1-4	0.4	8.8	0.9245	−3.052	Keratin associated protein 1-4
8172776	–	−0.4	5.4	0.9245	−3.057	ENSG00000197185
7944876	NRGN	−0.4	10.7	0.9245	−3.084	Neurogranin (protein kinase C substrate, RC3)
8112274	ELOVL7	−0.6	9.6	0.9245	−3.109	ELOVL family member 7, elongation of long chain fatty acids (yeast)
8041711	RNU6-2	0.4	4.2	0.9245	−3.117	RNA, U6 small nuclear 2
7961026	OVOS2	−1.0	6.5	0.9245	−3.143	Ovostatin 2
7905481	–	0.6	5.2	0.9245	−3.168	ENSG00000183586
8144719	LOC100271832	0.8	3.1	0.9245	−3.169	Hypothetical LOC100271832
8162529	–	0.4	6.3	0.9245	−3.178	ENSG00000209183
7982088	SNRPN	0.6	3.1	0.9245	−3.181	Small nucleolar RNA, C/D box 115-8
8130662	–	0.4	4.3	0.9245	−3.203	ENSG00000208409
7926889	LYZL1	0.5	4.0	0.9245	−3.205	Lysozyme-like 1
7913976	–	0.4	3.6	0.9245	−3.214	ENSG00000209703
7963208	–	0.4	3.5	0.9245	−3.222	ENSG00000210200
7915252	BMP8B	0.4	7.5	0.9245	−3.223	Bone morphogenetic protein 8b
8178205	HLA-DQA2	0.5	7.9	0.9245	−3.233	Major histocompatibility complex, class II, DQ alpha 1
8116313	LOC441120	0.4	6.1	0.9245	−3.249	Hypothetical protein LOC729569
7892934	–	−0.6	7.1	0.9245	−3.275
8168687	RN5S9	0.5	4.7	0.9245	−3.280	RNA, 5S ribosomal 9
8063716	TUBB1	−0.5	11.6	0.9245	−3.293	Tubulin, beta 1
8111772	DAB2	−0.5	8.7	0.9245	−3.297	Disabled homolog 2

“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. .

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) 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.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.