Exposure to workplace bullying, microRNAs and pain; evidence of a moderating effect of miR-30c rs928508 and miR-223 rs3848900

Abstract

Prolonged exposure to bullying behaviors may give rise to symptoms such as anxiety, depression and chronic pain. Earlier data suggest that these symptoms often are associated with stress-induced low-grade systemic inflammation. Here, using data from both animals and humans, we examined the moderating role of microRNAs (miRNAs, miRs) in this process. In the present study, a resident-intruder paradigm, blood samples, tissue harvesting and subsequent qPCR analyses were used to screen for stress-induced changes in circulating miRNAs in rats. The negative acts questionnaire (NAQ), TaqMan assays and a numeric rating scale (NRS) for pain intensity were then used to examine the associations among bullying behaviors, relevant miRNA polymorphisms and pain in a probability sample of 996 Norwegian employees. In rats, inhibited weight gain, reduced pituitary POMC expression, adrenal Nr3c1 mRNA downregulation, as well as increased miR-146a, miR-30c and miR-223 in plasma were observed following 1 week of repeated exposure to social stress. When following up the miRNA findings from the animal study in the human working population, a stronger relationship between NAQ and NRS scores was observed in subjects with the miR-30c GG genotype (rs928508) compared to other subjects. A stronger relationship between NAQ and NRS scores was also seen in men with the miR-223 G genotype (rs3848900) as compared to other men. Our findings show that social stress may induce many physiological changes including changed expression of miRNAs. We conclude that the miR-30c GG genotype in men and women, and the miR-223 G genotype in men, amplify the association between exposure to bullying behaviors and pain.

Lay summary

• Using an animal model of social stress, we identified miR-146a, miR-30c and miR-223 as potentially important gene regulatory molecules that may be involved in the stress response. Interestingly, human genotypes affecting the expression of mature miR-30c and miR-223 had a moderating effect on the association between exposure to bullying and pain. Subjects with the miR-30c rs928508 GG genotype had a significantly stronger association between exposure to bullying behaviors and pain than other subjects. The same was observed in men with the miR-223 rs3848900 G genotype, as compared to other men.

Keywords:

• Social stress
• miRNA
• Pain
• Rat
• Human
• Genotype
• SNP

Background

Exposure to negative social acts at the workplace by one’s peers or superiors is a common social stressor in contemporary working life with a global prevalence of about 15% (Agervold, 2007; Nielsen, Nielsen, Notelaers, & Einarsen, 2015). As previously described, the term “workplace bullying” refers to a situation in which a person repeatedly is subject to negative social acts and is unable to defend him/herself (Einarsen & Nielsen, 2015). Although there is no definitive list of bullying behaviors, bullying may involve actions like verbal harassment, spreading of rumors, physical confrontations, and social exclusion (Notelaers, Van der Heijden, Guenter, Nielsen, & Einarsen, 2018). Being a target of bullying at the workplace may have severe and detrimental consequences for well-being and workability.

A consistent body of evidence shows that bullying is associated with a subsequent increase in health complaints (Nielsen, Mageroy, Gjerstad, & Einarsen, 2014; Verkuil, Atasayi, & Molendijk, 2015), including self-reported pain (Jacobsen, Nielsen, Einarsen, & Gjerstad, 2018; Kaaria, Laaksonen, Rahkonen, Lahelma, & Leino-Arjas, 2012). While previous research has established the importance of cognitive factors as in the mechanisms underlying the association between bullying and somatic health outcomes (Mikkelsen & Einarsen, 2002; Nielsen, Matthiesen, & Einarsen, 2008), less is known about the role of biological factors. Especially, we think there is a need for more knowledge that enlighten the interaction between social stress and genetic factors with regard to pain.

Previous data show that social stress may be associated with activation of the anterior cingulate cortex (Eisenberger, Lieberman, & Williams, 2003) that plays a crucial part in central pain processing (Buffington, Hanlon, & McKeown, 2005; Zhao et al., 2006). Also, several lines of evidence show that persistent exposure to psychosocial stress leads to activation of microglia in other stress-responsive brain regions – including the prefrontal cortex, hypothalamus and amygdala (Hinwood et al., 2013; Johnson et al., 2005; Tynan et al., 2010; Wohleb et al., 2011). Therefore, exposure to bullying may influence pain signaling through central neuro-inflammatory processes (Dafny et al., 1996; Gracely, Petzke, Wolf, & Clauw, 2002; Manning & Mayer, 1995).

Social stress may also elicit changes in the autonomic nervous system (McQuade & Breaux, 2017; Sgoifo et al., 1999) and the hypothalamic–pituitary–adrenal gland (HPA) (Weaver et al., 2004). Previous data show that repeated social defeat may lead to increased myelopoiesis in the bone marrow and an accumulation of myeloid immune cells in the spleen (Engler, Bailey, Engler, & Sheridan, 2004). Moreover, an upregulation of glucocorticoid-insensitive circulating monocytes following social isolation has been reported (Cole et al., 2015). Thus, stress-induced low grade systemic inflammation may be associated with pain conditions like persistent low back pain (Pedersen, Schistad, Jacobsen, Roe, & Gjerstad, 2015), chronic widespread pain (Stensson, Ghafouri, Gerdle, & Ghafouri, 2017) and fibromyalgia (Backryd, Tanum, Lind, Larsson, & Gordh, 2017).

Recent data show that inflammatory processes may be modulated by miRNAs binding to the complementary 3′-UTR regions of mRNA molecules (Jonas & Izaurralde, 2015). For instance, miR-155 downregulates negative regulators of myeloid cell proliferation, resulting in increased immune cell numbers (O'Connell, Chaudhuri, Rao, & Baltimore, 2009), whereas mir-146a dampens the expression of pro-inflammatory mediators, such as IL-1β, TNF and IL-6 (Boldin et al., 2011). Interestingly, β-adrenergic and glucocorticoid signaling may affect the miRNA expression profile in target tissues (Hou et al., 2012; Smith, Shah, & Cidlowski, 2010), suggesting that the neuro-immune interface may involve regulatory miRNAs.

In the present study, an animal model was used to screen for stress-associated circulating miRNAs. The miRNA screening was then followed up in a human cohort, where we examined the potential moderating effect of single nucleotide polymorphisms (SNPs) affecting the expression of the most relevant miRNAs.

Methods

Animals

A resident-intruder paradigm, where intruders were exposed to a dominant resident one hour each day for one week, was used to study stress-induced changes in the expression of circulating miRNA. Each of the 10 male long Evans rats (500–550 g) used as residents were housed with a female long Evans rat (250 g). The 10 male Sprague–Dawley rats (400–500 g) used as intruders were housed in pairs, as were the 10 male Sprague–Dawley rats (400–500 g) used as controls (Janvier Labs, Le Genest St Isle – France). The different strains were kept in separate rooms. The experiments were performed in the dark period of an artificial 12 h light/12 h dark cycle. All animal experiments were approved by the Norwegian Food Safety Authority and performed in conformity with the laws and regulations controlling experiments and procedures on live animals in Norway.

The stress conditioning was performed by temporarily removing the female rat and introducing the intruder animal into the resident cage. The rats were separated immediately following violent physical contact. When separated, the cage was divided into two compartments by a perforated plastic wall, allowing the intruder to see, smell and hear the resident rat. After 60 min in the resident cage, the intruder rat was returned to its home cage. This procedure was repeated for 7 days, with the intruder animals being introduced to a new resident animal every day. The control animals followed the same procedure except that they visited a foreign cage without a resident rat. All animals were weighed 1 h following the conditioning.

Tissue harvesting

On the last day of the experiment, all Sprague–Dawley rats were euthanized by dislocation of the neck under isoflorane anesthesia. The pituitary and adrenal glands were harvested, frozen on liquid nitrogen and later stored in a −80 °C freezer.

Blood sampling

Blood samples were taken from intruder and control animals at baseline – 1 week before conditioning – and after conditioning at day 0, 3 and 7. In order to minimize stress during sampling, the animals were anesthetized with hypnorm (0.75 mL/kg, s.c., fentanyl citrate 0.315 mg/ml, flanisone 10 mg/ml, methyl parahydroxybenzoate 1.80 mg/ml and propyl parahydroxybenzoate 0.20 mg/ml). The vena saphena was then punctured with a small hypodermic needle and 500 µL blood was collected from each animal in a microvette tube coated with EDTA and spun down at 2000 rcf for 5 min. The plasma fraction and cell fraction were separately frozen on liquid nitrogen and stored in a −80 °C freezer.

RNA isolation and cDNA synthesis

The AllPrep DNA/RNA/miRNA Universal Kit (Qiagen) was used to isolate total RNA from the frozen pituitary and adrenal samples. Total RNA concentrations were measured using the Nanodrop 8000 system (Thermo Fischer). Synthesis of cDNA from these tissues was carried out using the qScript microRNA cDNA Synthesis Kit (Quanta) with 1 µg of total RNA as input. As previously described (Moen et al., 2017), total RNA was extracted from 100 µL plasma using the miRNeasy serum plasma isolation kit (Qiagen) following the manufacturer’s protocol. In order to adjust for variations in RNA extractions, synthetic C-elegans (C) mir-39-3p (Qiagen) was spiked in at a working concentration of 6 × 10⁷ copies/µL after the initial denaturation. Total RNA was eluted in 14 µl RNase-free water. A fixed volume of 7 µl eluate was used as input for the cDNA synthesis using the miScript II RT kit (Qiagen).

Tissue mRNA analysis

As previously described (Egeland, Moen, Pedersen, Brisby, & Gjerstad, 2013), the qPCR analysis was performed on a StepOnePlus qPCR machine (Applied Biosciences, USA). Primers were designed using Primer Express 2.0 (Applied Biosystems, Foster City, CA, USA) and checked for specificity by performing a BLAST search. Effort was made to design primers without nonspecific binding (the melting curves indicated no bi-products). For more details about the primers (Sigma Aldrich), see Supplementary Table 1. Target genes were normalized to β-actin as internal reference.

Plasma miRNA analysis

The Rat miRNome miScript miRNA PCR Array containing 653 miRNAs (cat.no. MIRN-216Z, Qiagen) was used to measure miRNA expression in plasma of the rats. All miRNAs with stable C_t values below 28 were examined in the screening. SNORD68 and miR-39-3p were used for normalization of miRNA expression. First, using pooled cDNA samples, the differences in miRNA levels from baseline to day 7 in both groups were calculated. The miRNAs were then listed according to difference in fold change (intruder group versus control group). From this list, the 5 miRNAs with the highest difference in fold change were selected for further analysis at all four time points, using all samples. All primers were delivered by Qiagen. qPCR was performed with 2x QuantiTect SYBR Green PCR Master Mix (Qiagen) on a StepOnePlus qPCR machine (Applied Biosciences, Foster City, CA, USA) with the following conditions: 95 °C for 15 min, followed by 40 cycles at 94 °C for 15 s, 55 °C for 30 s, and 70 °C for 30 s.

Human cohort

A random sample of 5000 employees from the Norwegian working force was drawn from The Norwegian Central Employee Register by Statistics Norway. As previously described (Jacobsen et al., 2018), sampling criteria were adults between 18 and 60 years of age employed in a Norwegian enterprise. Questionnaires were distributed through the Norwegian Postal Service during the spring 2015, with a response rate of 32%. Informed consent was given by the respondents and subjects who gave consent were also sent the saliva collection kits. The number of subjects included in the regression analyses ranged from 414 to 996. The survey was approved by the Regional Committee for Medical Research Ethics for Eastern Norway.

Instruments

Exposure to bullying behaviors in the workplace over the last 6 months was measured with the nine-item version of the Negative Acts Questionnaire – Revised (NAQ-R) inventory (Einarsen, Hoel, & Notelaers, 2009; Notelaers et al., 2018). NAQ-R describes negative and unwanted behaviors that may be perceived as bullying if occurring on a regular basis. The NAQ-R contained items referring to both direct (e.g. openly attacking the victim) and indirect (e.g. social isolation, slander) behaviors. It also contained items referring to personal as well as work related forms of bullying. For each item the respondents were asked how often they had been exposed to the behavior at their present worksite during the last six months. Response categories range from 1 to 5 (“never,” “now and then,” ”monthly,” “weekly” and “daily”). Cronbach’s alpha for the NAQ-R was 0.86.

To assess pain, subjects were asked to rate their mean pain intensity throughout the last week using an 11 point (0–10) numeric rating scale (NRS) with endpoints “no pain” and “worst imaginable pain.”

Genotyping

Collection of saliva and extraction of genomic DNA was done using OrageneRNA sample collection kit (DNA Genotech Inc. Kanata, Ontario, Canada) according to the manufacturer’s instructions. Quantification of genomic DNA was performed using the Nanodrop 8000 system (Thermo Fischer). SNP genotyping was carried out using custom TaqMan SNP genotyping assays (Applied Biosystems, Foster City, CA, USA). Approximately 45 ng genomic DNA was amplified in a 5 µl reaction mixture in a 384-well plate containing 2× TaqMan genotyping master mix (Applied Biosystems Foster City, CA, USA) and 1× assay mix, the latter containing the respective primers and probes. The probes were labeled with the reporter dye FAM or VIC to distinguish between the two alleles. Samples were amplified on a Quantstudio 5 machine (Applied Biosystems, Foster City, CA, USA) following an initial denaturing step at 95 °C for 10 min. The amplification consisted of 40 cycles of denaturing at 95 °C for 15 s and annealing and elongation at 60 °C for 60 s. Approximately 10% of the samples were re-genotyped and the concordance rate was 100%.

Statistical analyses

Differences in the change in bodyweight over time between intruder rats and control rats were analyzed using a two-way rmANOVA with a Greenhouse–Geisser correction. Pituitary and adrenal gland gene expression was investigated using fold change values for each sample, defined as the expression of the target gene normalized to the expression of β-actin and the mean level in the control group. Group differences were then compared using the Mann–Whitney rank sum test.

Data analysis on the miRNome array was performed using the Web-based software delivered by Qiagen. In brief, the ΔC_t-value for each miRNA profiled in the plate was calculated using the formula ΔC_t = C_t^miRNA – C_t^{SNORD68,mir-39-3p}, where C_t^{SNORD68,mir-39-3p} represented the mean of SNORD68 and mir-39-3p. For each miRNA, ΔΔC_t between the two time points for both groups was calculated using the formula: ΔΔC_t = ΔC_t (Day 7) −ΔC_t (Baseline). The fold-change was calculated as 2^−ΔΔCt.

The top 5 miRNAs based on difference in fold change between the control group and the stress group where followed up. Fold change values for each sample were defined by the expression of the target miRNA normalized to the expression of the reference genes SNORD68 and mir-39-3p. The abundance of each miRNA was calculated by the comparative C_t method and 2^−ΔCt values were then normalized to the mean value of the group at baseline. A false discovery rate (FDR) approach, with p values from each of the top 5 miRNAs plus assumed uniform distribution of p-values between 0.05 and 1.00 for other miRNAs, was used to correct for multiple testing. The top 3 miRNAs, ranked by q-value, were then followed up in the human cohort.

The mean-score of the nine items in the NAQ-R was used to measure exposure to bullying. To explore the hypotheses about main and moderating effects, we conducted hierarchical regression analyses to test for linear associations between exposure to bullying behaviors and pain intensity, as well as the interactive effects of exposure to bullying with the different genotypes, with regard to the same outcome. In order to examine the moderating role of the genotypes, we followed the recommendations for interaction analyses provided by Baron and Kenny (Baron & Kenny, 1986). Genotypes were included as categorical moderators using one of the homozygous genotypes as reference group. For SNPs where the number rare-rare homozygote subjects were lower than 50, these were grouped together with the rare-common heterozygote subjects. The interaction analysis was conducted in two steps. Control variables, exposure to bullying and the genotypes were entered as predictors in the first step, whereas the interaction term (exposure to bullying*genotype) was entered in the second step. A significant interaction term in the second step was considered as an interaction effect.

As the scores on the NAQ (skewness: 3.92; kurtosis: 22.02) were non-normally distributed, all analyses were conducted using bootstrapping (5000 resamples). Statistical analyses were conducted with Stata 14 (StataCorp). The level of significance was set to p < .05.

Results

Animal study – screening of miRNAs

The resident-intruder paradigm caused a significant attenuation of weight gain in the stress-exposed rats (Figure 1(A)). A significant decrease in the expression of proopiomelanocortin (POMC, 0.49 ± 0.13 fold of control), but no significant change in the expression of Nr3c1 (Figure 1(B,C)) were demonstrated in the pituitary gland. We also observed decreased expression of Mcr2 and Nr3c1 in the adrenal gland, but only the Nr3c1 decrease reached significance (0.33 ± 0.03 fold of control, Figure 1(D,E)). Regarding the miRNA screening, the top 3 miRNAs, ranked by q-value, were miR-146a, miR-30c and miR-223 (Table 1, Figure 2(A–C)).

Figure 1. (A) Bodyweight of exposed animals compared to controls at baseline, day 0, day 3 and day 7 of the stress paradigm, two-way rmANOVA, df = 1.209, p = .025. (B–E) Gene expression of the ACTH precursor POMC and the glucocorticoid receptor Nr3c1 in the pituitary gland, and Mcr2 and Nr3c1 in the adrenal gland following 7 days of social stress. Mann–Whitney rank sum test, *p < .05, **p < .01. All data are given as the mean ± SEM.

Figure 2. Fold expression in plasma after 0, 3 and 7 days of social stress exposure of (A) miR-146a, (B) miR-30c and (C) miR-223. Mann–Whitney rank sum test at day 7, using a false discovery rate (FDR) approach, q = 0.024, q = 0.048 and q = 0.144 for (A), (B) and (C), respectively. Data are given as the mean ± SEM.

Table 1. List of top 5 miRNAs from the screening (animal model), ranked by q-value; the top 3 miRNAs were followed up in the human study (working population).

MicroRNA (miRNA) ID	Fold change exposure	Fold change control	Difference in fold change	Mann–Whitney rank sum test p-value	False discovery rate (FDR) q-value
rno-miR-146a-5p	2.148	−0.334	2.489	.001	0.024
rno-miR-30c-5p	1.761	−0.437	2.199	.004	0.048
rno-miR-223-3p	1.019	−0.494	1.514	.016	0.127
rno-miR-320-3p	1.888	−0.098	1.986	.256	0.651
rno-miR-126a-3p	1.857	−0.033	1.891	.307	0.672

Human study – role of miR-146a, miR-30c and miR-223 in subjects exposed to bullying

As expected, a clear positive correlation between the NAQ and NRS scores was seen. However, no interaction between the NAQ score and the miR-146a rs2910164 genotype with regard to NRS was demonstrated (Supplementary Table 2, Table 2 and Figure 3(A)). Still, the data revealed a significantly stronger relationship between the NAQ and NRS scores in subjects with the miR-30c GG rs928508 genotype compared to subjects with GA and AA genotype. A clear allele dependent effect was observed: individuals with GG reported more pain following bullying than individuals with GA, which in turn reported more pain than AA (Supplementary Table 3, Table 3 and Figure 3(B)). Inclusion of the interaction term increased the explained variance from 8.1% to 9.0%. Since miR-223 is located on the X-chromosome, the data of men and women were analyzed separately. A significantly stronger relationship between the NAQ and NRS scores was demonstrated in men with the miR-223 rs3848900 G genotype compared to men with the A genotype (Supplementary Table 4, Table 4 and Figure 3(C)). Inclusion of the interaction term also increased the explained variance from 11.3% to 12.4%. No such relationship was observed in women (Supplementary Tables 5 and 6, Supplementary Figure 1).

Figure 3. The relationship between NAQ and pain intensity (NRS). Subjects were divided into groups based on (A) miR-146a genotype: GG and GC/CC (used as reference), (B) miR-30c genotype: AA, AG and GG (used as reference) and (C) miR-223 genotype among men: A and G (used as reference). Age, tobacco use and education were included as covariates in all analyses, and sex was included as covariate in (A) and (B).

Table 2. Hierarchical regression with miR-146a genotype GC/CC as reference (bootstrapping with 5000 resamples).

	Pain	B	Std. Err	p Value	95% conf. interval
Step 1	Rsq: 0.0817
	Age	0.008	0.007	.243	−0.005 to 0.021
	Sex	0.569	0.140	.000	0.295 to 0.843
	Tobacco	0.454	0.180	.012	0.102 to 0.806
	Education
	High school	−0.221	0.298	.459	−0.804 to 0.363
	University <4y	−0.807	0.286	.005	−1.367 to −0.247
	University >4y	−1.220	0.291	.000	−1.790 to −0.651
	miR-146a
	GG	0.226	0.142	.111	−0.046 to 0.499
	NAQ	1.002	0.259	.000	0.494 to 1.510
Step 2	Rsq: 0.0818
	Age	0.008	0.007	.245	−0.005 to 0.022
	Sex	0.569	0.140	.000	0.295 to 0.842
	Tobacco	0.456	0.180	.011	0.104 to 0.808
	Education
	High school	−0.217	0.297	.465	−0.800 to 0.366
	University <4y	−0.806	0.286	.005	−1.366 to −0.246
	University >4y	−1.220	0.290	.000	−1.789 to −0.651
	miR-146a
	GG	0.076	0.646	.904	−1.165 to 1.318
	NAQ	0.926	0.437	.034	0.069 to 1.782
	miR-146a x NAQ
	GG	0.126	0.535	.814	−0.923 to 1.174

Table 3. Hierarchical regression with genotype GG as reference (bootstrapping with 5000 resamples).

	Pain	B	Std. Err	p Value	95% conf. interval
Step 1	Rsq: 0.081
	Age	0.011	0.007	.109	−0.003 to 0.025
	Sex	0.559	0.144	.000	0.277 to 0.840
	Tobacco	0.395	0.186	.034	0.030 to 0.760
	Education
	High school	−0.204	0.309	.509	−0.810 to 0.402
	University <4y	−0.822	0.294	.005	−1.398 to − 0.246
	University >4y	−1.205	0.297	.000	−1.787 to −0.623
	miR-30c
	AG	−0.222	0.210	.290	−0.633 to 0.189
	AA	−0.286	0.219	.191	−0.714 to 0.143
	NAQ	0.987	0.266	.000	0.466 to 1.509
Step 2	Rsq: 0.090
	Age	0.011	0.007	.112	−0.003 to 0.025
	Sex	0.561	0.143	.000	0.281 to 0.841
	Tobacco	0.420	0.186	.024	0.056 to 0.784
	Education
	High school	−0.160	0.303	.598	−0.754 to 0.435
	University <4y	−0.773	0.289	.008	−1.340 to −0.206
	University >4y	−1.183	0.292	.000	−1.755 to −0.611
	miR-30c
	AG	1.846	0.972	.057	−0.059 to 3.750
	AA	2.649	0.995	.008	0.698 to 4.600
	NAQ	2.930	0.772	.000	1.417 to 4.444
	DNMT3B x NAQ
	AG	−1.791	0.851	.035	−3.459 to −0.122
	AA	−2.507	0.864	.004	−4.200 to −0.815

Table 4. Hierarchical regression of men with miR-223 genotype G as reference (bootstrapping with 5000 resamples).

	Pain	B	Std. Err	p Value	95% conf. interval
Step 1	Rsq: 0.113
	Age	0.026	0.010	.007	0.0087 to 0.046
	Tobacco	0.620	0.216	.004	0.197 to 1.043
	Education
	High school	−0.196	0.395	.620	−0.970 to 0.578
	University <4y	−1.049	0.381	.006	−1.796 to −0.302
	University >4y	−1.093	0.390	.005	−1.858 to −0.328
	miR-223
	A	−0.017	0.255	.946	−0.518 to 0.484
	NAQ	1.314	0.398	.001	0.535 to 2.093
Step 2	Rsq: 0.124
	Age	0.027	0.010	.005	0.008 to 0.046
	Tobacco	0.603	0.218	.006	0.175 –1.030
	Education
	High school	−0.144	0.389	.711	−0.907 to 0.619
	University <4y	−1.041	0.378	.006	−1.781 to −0.300
	University >4y	−1.070	0.388	.006	−1.829 to −0.310
	miR-223
	A	2.167	1.097	.048	0.016 to 4.318
	NAQ	2.774	0.795	.000	1.125 to 4.332
	miR-146a x NAQ
	A	−1.813	0.916	.048	−3.607 to 0.018

Discussion

Previous studies have shown that prolonged exposure to social stress in mammals may induce behavioral (Rygula et al., 2005; Zelena, Haller, Halasz, & Makara, 1999) and/or metabolic (Liu et al., 2016) changes. Therefore, such stress probably affects both food intake and energy requirements in social animals. In accordance with these findings, our data showed that the resident-intruder paradigm caused reduced weight gain in the intruder animals, demonstrating severe detrimental consequences of social stress.

Earlier findings show that prolonged stress may reduce the expression of the ACTH precursor POMC in the pituitary gland (Chen, Tang, & Yang, 2008; Hauger, Millan, Lorang, Harwood, & Aguilera, 1988). Several signaling pathways may be responsible for these changes. For instance, vasopressin signaling from the hypothalamus (Bartanusz et al., 1993; Levin, Blum, & Roberts, 1989), and elevated circulating glucocorticoid levels (Young, Kwak, & Kottak, 1995) may downregulate POMC. In line with these findings, the present data showed that 1 week of repeated social stress caused reduced POMC expression in the pituitary gland. However, upregulation of POMC following prolonged stress exposure has also been reported (Kiss & Aguilera, 1993; Lopez-Calderon, Ariznavarreta, & Chen, 1991).

Like other components of the HPA axis, the adrenal gland expresses glucocorticoid receptors (Ceccatelli et al., 1989). The hypothalamic paraventricular nucleus (Evanson, Tasker, Hill, Hillard, & Herman, 2010; Tasker & Herman, 2011) and the pituitary gland (Deng et al., 2015) may be major players of the HPA axis negative feedback mechanism. In addition, the adrenal gland may also have a local negative feedback loop that controls glucocorticoid release (Riester et al., 2012). Interestingly, glucocorticoid signaling in the adrenal gland may participate in the regulation of adrenal androgens (Paust et al., 2006), which in turn regulates glucocorticoid production (Pinto et al., 2015). Thus, the observed downregulation of glucocorticoid receptors in the adrenal gland suggests weakened feedback inhibition of glucocorticoid production, potentially leading to prolonged and enhanced responses following further stressful stimuli.

Many miRNAs, including those emphasized in the present study, may act as anti-inflammatory regulators. Notably, miR-146a and miR-223 are both enriched in glial cells (Jovicic et al., 2013), where they negatively regulate toll-like receptor signaling (Taganov, Boldin, Chang, & Baltimore, 2006; Wang et al., 2015). Moreover, miR-146a may reduce IL-1 dependent inflammation (Gu et al., 2015) and has been shown to decrease the expression of several other pro-inflammatory mediators, such as TNF, IL-6 (Nahid, Pauley, Satoh, & Chan, 2009) and IL-8 (Bhaumik et al., 2009).

Evidence exists that miR-146a may be downregulated in the prefrontal cortex of patients with depression (Smalheiser et al., 2012) and in the monocytes of patients with bipolar disorders (Weigelt et al., 2013). However, in our human cohort, the miR-146a rs2910164 SNP – located 60 basepairs downstream of the first nucleotide of pre-miR-146a, and which may influence the expression of mature miR-146 (Jazdzewski et al., 2008) – did not influence vulnerability toward bullying behaviors with regard to pain.

Downregulation of miR-30c may induce systemic inflammation by regulating pro-atherosclerosis pathways in macrophages (Ceolotto et al., 2017). Moreover, it has been reported that miR-30c may exert important epigenetic action by selectively downregulating DNA methyltransferase 3B (DNMT3B) (Liu et al., 2012). In general, DNA methylation (CH3 groups added to cytosine bases in CpG islands of the DNA) is associated with gene repression. Interestingly, miR-30c rs928508 is located in position bp +419, not very far from the coding region of the mature miR. Previous data show that miR-30c rs928508 G in this “flanking region” of miR-30c may inhibit the transition from pri-miRNA to pre-miRNA, ultimately leading to lower levels of mature miR-30c (Chen, Liu, Hu, & Shen, 2012; Hu et al., 2011). It is therefore likely that this genetic variant reduces the miR-30c expression and upregulates DNMT3B.

Previous data show that the miR-30 family may affect cell adhesion, cell cycle, stress response and EGF activation (Izzotti et al., 2009). Here, we show a clear interaction between exposure to negative behavior and miR-30c genotype; this A>G polymorphism clearly influenced vulnerability to bullying behaviors, i.e. pain sensitivity. One explanation of the present finding that GG subjects report more pain following bullying than GA subjects, that in turn report more pain than the AA subjects, could be that the miR-30c A allele has a protective effect. Thus, genetic variability in the gene encoding miR-30c may affect the response to bullying behaviors.

In the periphery, miR-223 seems to be involved in regulating the proliferation and activity of myeloid cells, including neutrophils (Bauernfeind et al., 2012; He et al., 2017; Johnnidis et al., 2008) and macrophages (Chen et al., 2012). Moreover, miR-223 released in exosome-like vesicles from nucleus pulposus cells may regulate inflammation following disc herniation (Moen et al., 2017). A protective role of miR-223 against inflammation-driven pain states, such as atherosclerosis (Wang et al., 2015) and chronic lumbar radicular pain (Moen et al., 2017), has also previously been suggested.

In addition, animal experiments show that miR-223 could influence central pain processing by targeting AMPAR and NMDAR subunits, thereby reducing neuronal excitability in response to glutamate (Harraz, Eacker, Wang, Dawson, & Dawson, 2012). Recent evidence shows that the expression of miR-223 may be regulated through a C/EBPa binding site between position bp-730 and -709 (Fazi et al., 2005). In our human cohort, we demonstrate that the miR-223 rs3848900 G allele within this binding site was associated with increased bullying-induce pain in men. The fact that miR-223 is located on the x-chromosome explains why we see this association in men, but not in women.

Conclusions

Our findings support the hypothesis that stress-induced neuronal or inflammatory processes may be associated with pain. Subjective health complaints in individuals exposed to bullying behavior was moderated by the miR-30c and miR-223 genotype. Interestingly, the miR-30c G allele in general, and the miR-223 G allele in men, appeared to increase the risk of pain in subjects exposed to workplace bullying. Thus, the present data show that genetic factors are important for the mechanisms underlying the association between bullying and somatic health outcomes.

Moreover, the present data show that the association between exposure to negative social acts among vulnerable individuals might be more potent than previously assumed. Employers, organizations, and health professionals need to acknowledge that individual biological differences should be taken into consideration when developing interventions against bullying at the workplace. We conclude that stress-induced neuronal or inflammatory processes influenced by regulatory miRNAs may be important to our future understanding of emotional and physical reactions following prolonged exposure to bullying behaviors.

Ethics approval

All animal experiments were approved by the Norwegian Food Safety Authority (application ID: 8212). The human arm of the study was approved by the Norwegian Regional Committee for Medical Research Ethics (REK 2014/1725)

Notes on contributions

DPJ, MBE, and JG; AM, DPJ, MBE, DR and JG performed the research; DPJ, MBE, IN, DR, MBN, SE and JG analyzed the data; DPJ and JG wrote the paper. All authors read and approved the final manuscript.

Acknowledgements

The authors thank Fang-Chin Lin, Aqsa Mahmood, Anne-Mari Gjestvang Moe, Tiril Schjølberg and Øivind Skare for excellent work.

Availability of data and materials

All data generated or analyzed during this study are included in this published article (and its supplementary information files).

Disclosure statement

No potential conflict of interest was reported by the authors.

Correction Statement

This article has been republished with minor changes. These changes do not impact the academic content of the article.

Additional information

Funding

This work was funded by the Norwegian research council [grant number 237777/250127].