Diosgenin Glucoside Protects against Spinal Cord Injury by Regulating Autophagy and Alleviating Apoptosis

Spinal cord injury (SCI) is a severe traumatic lesion of central nervous system (CNS) with only a limited number of restorative therapeutic options. Diosgenin glucoside (DG), a major bioactive ingredient of Trillium tschonoskii Max., possesses neuroprotective effects through its antioxidant and anti-apoptotic functions. In this study, we investigated the therapeutic benefit and underlying mechanisms of DG treatment in SCI. We found that in Sprague-Dawley rats with traumatic SCI, the expressions of autophagy marker Light Chain 3 (LC3) and Beclin1 were decreased with concomitant accumulation of autophagy substrate protein p62 and ubiquitinated proteins, indicating an impaired autophagic activity. DG treatment, however, significantly attenuated p62 expression and upregulated the Rheb/mTOR signaling pathway (evidenced as Ras homolog enriched in brain) due to the downregulation of miR-155-3p. We also observed significantly less tissue injury and edema in the DG-treated group, leading to appreciable functional recovery compared to that of the control group. Overall, the observed neuroprotection afforded by DG treatment warrants further investigation on its therapeutic potential in SCI.

Introduction
Spinal cord injury (SCI) remains a rare disease in comparison to traumatic brain injury, concussion, or stroke in the younger population, but it can cause permanent disability or loss of movement and sensation. It can lead to a series of complications with complex pathological process, and more difficulties in rehabilitation, with only a limited number of therapeutic options to improve functional recovery. In China, the incidence of SCI is increasing yearly [1]. Therefore, SCI repair is a hot topic in clinical research.
At present, the pathological and repairing mechanisms of SCI mainly focuses on the prevention and reversible regulation of secondary injury of spinal cord, which will provide reliable clinical treatment strategies for the recovery and regeneration of neurons after SCI. The initial trauma is followed by complex, poorly understood changes resulting in further destruction of the spinal cord tissue over a period of several months, and consequently loss of the neural function. The apoptosis of neurons is an important mechanism in the pathogenesis of SCI [2,3]; therefore, the major strategy is

DG Decreased the Structural Damage of Spinal Cord Tissue and Promoted Functional Recovery after SCI
To evaluate the therapeutic effect of DG on locomotor recovery after SCI, functional recovery was evaluated after injury for 21 days using Basso, Beattie, and Bresnahan (BBB) scores. BBB scores of the sham group were about 21 normally, whereas the BBB scores of SCI group and DG-treated group were below the normal score after injury (p < 0.01) (Figure 1). There was no significant difference in BBB scores between SCI group and DG-treated group after SCI for four days; however, BBB scores revealed an increase after the contusion with DG treatment from the fifth day when compared with SCI model group (p < 0.05). These data indicate that DG may influence the functional improvement of locomotor activity after SCI. In general, functional recovery gradually improved both in diosgenin glucoside (DG) and SCI groups from the fourth day after SCI. Except for the first four days the experimental period, the DG group showed significantly improved BBB scores in comparison with the SCI group. ** p < 0.01 in the DG group when compared with the sham group. # p < 0.05 in the DG group when compared with the SCI group (n = 10 per group).
The difference in morphology of tissues from sham, SCI, and DG groups was analyzed by Hematoxylin-Eosin (HE) staining after the injury for one, seven, 14, and 21 days (Figure 2A). The ratio of the cavity area in the spinal cord cross-sectional area by HE staining was calculated somewhat, to histologically show the function of the spinal cord [28,29]. Severe damage of the dorsal white matter and central grey matter was observed in SCI rats when compared with the sham group. In contrast, obvious therapeutic effects were observed in DG-treated group, with a decreased cavity in the dorsal white matter and central gray matter, as compared with SCI group (p < 0.01) ( Figure 2D). (D,E) graphic presentation indicates the ratio of cavity area (cavity area, spinal cord cross-sectional area), the number of surviving neurons and the number of Nissl bodies. All data were expressed as M ± SD. ** p < 0.01 in the DG group when compared with the sham group. ## p < 0.01 in DG group when compared with the SCI group (n = 5 per group).
Neurons stained with Nissl were counted to examine the effect of DG on the loss of motor neurons in the ventral horn. The number of Nissl bodies was increased in the DG group when compared with that of the SCI group (p < 0.01) ( Figure 2B,E). In order to determine the morphological changes of the spinal cord and the structures of spinal cord from a sham operation, the SCI model, and DG treatment groups were examined by transmission electron microscopy (TEM). In the normal group, myelin sheaths were arranged in the integral structures and clear mitochondria. In comparison, the spinal cord in the SCI group revealed irregular-shaped, disordered structures with loose, twisted, wrinkled, or partially ruptured plates. Swollen and ruptured mitochondria were also observed in this group, with a greatly decreased number of mitochondria. However, in the DG group, myelin sheaths exhibited more regular shapes and integral structures, more uniformed lamellar structures, reduced degrees of loosening, improved degrees of rupture, and alleviated mitochondrial swelling, when compared with the SCI group ( Figure 2C).

DG Attenuated Apoptosis Caused by SCI
In order to evaluate DG-modulated cellular apoptosis, the expression levels of cleaved caspase-3, Bax, and Bcl-2 were analyzed by Western blot and immunohistochemistry (IHC). Western blot analysis of cleaved caspase-3 showed that DG significantly reversed the increased level of cleaved caspase-3 caused by SCI ( Figure 3B,E). Similarly, the expression of pro-apoptotic protein Bax was down-regulated, while the anti-apoptotic protein Bcl-2 was up-regulated after treatment with DG ( Figure 3C,E). Moreover, terminal-deoxynucleotidyl transferase (TdT)-mediated nick end labeling (TUNEL) staining of the spinal cord was performed to estimate the amount of apoptotic cells within in each group. As shown in Figure 3A,D, the number of apoptotic cells was increased significantly in the SCI group when compared with that of the sham group. In contrast, the number of TUNEL-positive neurons in the DG group was less than that of the SCI group (p < 0.01). Taken together, these results demonstrate that DG administration is an effective therapeutic strategy for the inhibition of apoptosis after SCI. (C,E) representative Western blots and quantification data of Bax, BCL-2, caspase-3, and β-actin in each group. Columns represent M ± SD for all data. * p < 0.05 and ** p < 0.01 in the DG group when compared with the sham group. # p < 0.05 and ## p < 0.01 in DG group when compared with the SCI group (n = 5 per group).

DG Down-Regulated miR-155-3p in Spinal Cord Tissue after SCI
By using a bioinformatics analysis, we examined two subtypes of miR-155, and the subtype miR-155-3p had altered expression after SCI. The expression of miR-155-3p was further confirmed by quantitative real-time polymerase chain reaction (qRT-PCR) analysis using the comparative Ct (∆∆Ct) method. The expression of miR-155-3p in spinal cord tissue was significantly different between groups. As shown in Figure 4, the expression of miR-155-3p was significantly increased in the SCI model group when compared with the sham group (p < 0.05). However, DG administration markedly mitigated the up-regulation of miR-155-3p after SCI (p < 0.05), suggesting that the up-regulation of miR-155-3p is associated with the SCI of rats. . The miR-155-3p was activated in spinal cord tissue from SCI rats. The expression of miR-155-3p in different groups was analyzed through quantitative real-time polymerase chain reaction (qRT-PCR). All data were expressed as M ± SD. ** p < 0.01 in the DG group when compared with the sham group. # p < 0.05 in the DG group when compared with the SCI group (n = 5 per group).

DG-Activated Rheb/mTOR Signal Pathway after SCI
To explore whether the Rheb/mTOR signal pathway was activated by DG, the expression of Rheb, mTOR, p-mTOR, AKT, p-AKT, and p70S6K was detected in spinal cord tissues from each group by Western blot analysis after SCI for seven days ( Figure 5A,B). Representative Western blots and quantification data of pAKT, AKT, p-mTOR, mTOR, p70S6K, Rheb, and β-actin in each group. All data were expressed as M ± SD. * p < 0.05 and ** p < 0.01 in the DG group when compared with the sham group. # p < 0.05 and ## p < 0.01 in the DG group when compared with the SCI group (n = 5 per group).
Rheb was significantly decreased, and mTOR and p70S6K was significantly reduced in spinal cord tissues from the SCI model group when compared with those in the sham group. DG treatment significantly up-regulated Rheb and p70S6K, and increased p-mTOR/mTOR ratio when compared with those in the SCI group. This result suggests that DG further activated the Rheb/mTOR signal pathway after SCI.

DG Rescued Abnormal Autophagy in Spinal Cord Tissues of SCI Rats
In order to verify the regulatory role of DG in autophagy, the expression of autophagy-associated proteins such as the ratio of LC3-II to LC3-I, Beclin1, and p62 was evaluated by Western blot and IHC. SQSTM1/p62, a polyubiquitin binding protein that is selectively incorporated into autophagosomes through direct binding to LC3 and efficient degradation during autophagy, was also examined. Thus, the total cellular levels of SQSTM1 reflect the extent of autophagy. The results showed that the ratio of LC3-II to LC3-I and the expression of Beclin1 were significantly decreased, and the expression of p62 was significantly increased in spinal cord tissue of SCI rats when compared with those in the sham group. However, the ratio of LC3-II to LC3-I and the expression of Beclin1 increased, while the expression of p62 decreased in DG group when compared with the SCI group. This suggests that dysfunctional autophagy in spinal cord tissue of SCI rats can be rescued by DG treatment ( Figure 6A-C). To further analyze the degree of autophagy in the SCI and DG groups at each designated time point (one, seven, 14, or 21 days), the expression of p62 and LC3 was evaluated by western blot. The expression levels of p62 and LC3 were significantly increased in a time-dependent manner at one, seven and 14 days, respectively, but decreased at 21 days in DG group when compared with the SCI group ( Figure 7A-D). These data suggest that autophagy is a dynamic process. DG treatment resulted in the reduction of p62 as compared with the SCI group, suggesting that DG increases autophagy by promoting degradation by autophagy. DG rescued the impaired autophagy in spinal cord tissues of rats with SCI. (A) Immunofluorescence staining of LC3 (red), GFAP (green) and DAPI-labeled nuclei (blue). Scale bars are 10 µm; (B,C) representative western blots and quantification data of LC3, Beclin1, p62 and β-actin in each group. All data were expressed as M ± SD. * p < 0.05 and ** p < 0.01 in the DG group when compared with the sham group. # p < 0.05 and ## p < 0.01 in the DG group when compared with the SCI group (n = 5 per group).

Discussion
SCI is a more important type of injury in clinical practice, and it is an important topic of clinical research in the field of medicine. At present, the pathological and repairing mechanisms of SCI mainly focus on the prevention and reversible regulation of secondary injury of spinal cord, which will prevent further death of neurons and glial cells in spinal cord tissue, and will provide reliable clinical treatment strategies for the repairing and regeneration of neurons after spinal cord injury. The treatment using bioactive ingredients from traditional Chinese herbs with low toxicity and side effects has become a novel research trend. Previous studies have confirmed that DG has neuroprotective effects against SCI. According to reports in the literature, the apoptosis of neurons, oligodendrocytes, astrocytes, and microglia cells in tissues from injured segments and adjacent segments from 4 hr to 3 weeks could be observed after SCI [30][31][32]. The apoptosis of spinal cord neurons caused the impairment and loss of spinal cord sensory and motor functions [29]. Consistent with our findings, it is known that trauma and ischemia of the spinal cord can cause neuronal apoptosis, which could be attenuated by DG treatment in the present study. Our data also indicates that DG may influence the functional improvement of locomotor activity after SCI, as confirmed by Basso, Beattie, and Bresnahan (BBB) scores and TUNEL staining, which could be an effective therapeutic strategy for inhibiting apoptosis after SCI. In the spinal cord of patients with SCI, apoptosis was mainly found in the lesion and its surrounding tissues, while the Waller change was observed in the lesion center. Similarly, DG treatment can reduce loosening degree, improve rupture degree, and alleviate mitochondrial swelling, as observed by TEM when compared with the SCI group.
Previous studies have shown that microRNAs can play a broad range of roles in cell differentiation, proliferation, and apoptosis by regulating the expression of hundreds of genes via RNA-induced silencing complexes (RISC) and targeting RNA expression processes [33]. MicroRNAs are rich in the central nervous system, and reveal a close relationship with the development of neurological diseases, as well as widely influence the signaling networks, leading to pathological responses after SCI. MicroRNAs have gained extensive attention, and could become a novel therapeutic target for the treatment of SCI. A large number of studies have shown that miR-155 is a broad regulator of inflammatory mediators, and plays a key regulatory role in immune and inflammatory disorders [34][35][36][37]. The development of neurodegenerative diseases is closely correlated with immunity and inflammation, suggesting that miR-155 may be a novel target for its treatment. A report [15] has documented that there is high expression with more than 2-fold change of miR-155 and other microRNAs (miRNAs) in brain tissue with multiple lesion locations. We previously demonstrated that DG treatment could downregulate the expression of miR-155 in a rat model of D-galactose-induced brain senescence [30]. Based on these studies and bioinformatics analysis, we analyzed the expression of miR-155 by qRT-PCR, demonstrating that miR-155-3p was aberrantly expressed after SCI. In addition, our findings have implicated that DG treatment could down-regulate miR-155-3p expression after SCI, suggesting that miR-155-3p may play an important role in the pathogenesis of SCI. Therefore, it is necessary to verify the target gene of miR-155-3p during DG treatment after SCI in order to further explore the underlying mechanisms of DG.
Rheb is regulated by the PI3K signaling pathway, including PKB and TSC1/2. Rheb modulates the mTOR complex and controls protein synthesis via the phosphorylation of 4E-BP1 and S6. Wan et al. [38] had demonstrated that Rheb, RICTOR, and RPS6KB2 (p70S6K) of the mTOR signal pathway are direct targets of miR-155. miR-155 can bind to the 3 -UTR of the target genes Rheb and p70S6K, thereby negatively regulating autophagy. Our results demonstrated that DG treatment activated the Rheb/p70S6K signal pathway. In addition, Rheb is beneficial for the attenuation of neurodegenerative diseases including Alzheimer's disease, Parkinson's disease, and Huntington disease by secreting brain-derived growth factors (BDGF) and glial-derived neurotrophic factors, suggesting that Rheb plays an important role in the central nervous system. The overexpression of Rheb in the nervous system can promote the growth of axons [39] and the production of myelin after birth, which plays a protective role in neurons [40][41][42]. Rheb is the upstream member of mTORC1, and plays an important role in the regulation of mTORC1 activity [43]. Rheb can bind to GTP to promote cell growth and cell cycle regulation, and insulin signaling pathways play an important role in maintaining the growth of organisms. mTOR kinase plays an important role in the transduction of extracellular to intracellular signals, and is an important target for the treatment of neurological diseases. The function of mTOR can be accomplished through a complex of mTORC1 and mTORC2. Rheb regulates the self-renewal of neural stem cells by regulating the activity of mTORC1, which plays an important role in the development of the cerebral cortex. In contrast, in injured neurons, Rheb is mainly expressed in astrocytes and neurons [44,45], indicating that Rheb is highly correlated with the proliferation and apoptosis of neural progenitor cells. Our data showed that DG treatment could significantly increase p-mTOR/mTOR ratio, implying that DG could further activate Rheb/mTOR signal pathway and inhibit miR-155-3p expression after SCI. In order to confirm whether miR-155-3p is related to the expression of Rheb, we will detect the expression of Rheb/mTOR signal pathway using a miR-155-3p mimic and inhibitor in vitro.
Autophagy is an important downstream result produced by the inhibition of the mTOR signal pathway, which can be used as a new therapeutic modality for SCI as previous report [46][47][48]. The activation of mTOR showed that p-mTOR/mTOR ratio in SCI rats was significantly increased when compared with the rats from the sham group, while DG treatment could significantly inhibit mTOR activation in SCI rats. Based on these results, it is plausible that DG can activate autophagy by activating mTORC1, which may be due to the impaired flux of autophagy during SCI progression. DG treatment could rescue dysfunctional autophagy and promote autophagy after SCI in a time-dependent manner, as illustrated by analyzing the degree of autophagy at each time point (one, seven, 14, or 21 days) through p62 and LC3 expression. Rheb expression was increased in light-induced retinal neurons, and the expression of caspase-3 was enhanced, with the overlap of Rheb and caspase-3 under the detection by double-labeled fluorescent immunoassay [49]. Rheb activation not only promotes cellular proliferation and differentiation, but also enhances cellular apoptosis in response to diverse toxic stimuli. Furthermore, we also found that the expression profiles of p62 and LC3 were parallel to Rheb in a time-dependent manner. DG treatment could enhance the expression of Rheb and decrease the expression of miR-155-3p, following the significant reduction of caspase-3 and Bax expression, and up-regulation of Bcl-2 expression. However, the potential mechanisms for Rheb activation in the presence of DG still need to be explored further. Therefore it is necessary to further clarify the molecular mechanisms and the potential relationship between the inhibition of SCI-induced apoptosis, and the activation of autophagy or the rescuing of impaired flux of autophagy during DG treatment after SCI, using corresponding autophagy activators and inhibitors in vitro, as well as to examine whether miR-155-3p/Rheb is relevant to neuronal apoptosis after SCI.
The outcomes of these studies will provide a reference for uncovering the partial underlying mechanisms and exploring the therapeutic targets of SCI and other autophagy-associated diseases during intervention, using natural products.

Animal Grouping and Treatments
A total of 120 healthy SD (Sprague Dawley) rats of SPF (Specific Pathogen Free) grade with body weights of 250 to 300 g were purchased from the Experimental Animal Center of Hubei province and housed in a temperature-controlled environment at room temperature of 20 to 24 • C and a relative humidity of 50 ± 20%. The rats were provided with standard feed and drinking water freely. All animal experiment protocols were reviewed and approved by the animal ethics committee of Hubei University for Nationalities (No. 2016(021)) in July 2016, and experimental animal welfare and ethical principles were followed.
A rat model with acute SCI was established according to the previously reported Allen's weight-striking method [50] with minor modification. After inducing anesthesia using 10% chloral hydrate solution at a dose of 3 mL/kg, the rats were fixed on the operation table in a prone position, and sequentially subjected to hair removal and disinfection. The location of T10 was determined, and an incision with a length of 2 to 3 cm was made at the center of T10. After cutting the tissues layer by layer, the processus spinosus of T9-T11 was fully exposed. The processus spinosus and lamina of T9-T11 were removed by clamp to expose the spinal canal and spinal dura. Using a modified spinal cord striking device, a balance weight of 5 g was made to fall freely at a height of 5 cm to hit the exposed spinal dura, thus resulting in SCI of T9-T11. The occurrence of tail-spasm swing and strong lower limb contraction was used as the indicator for evaluating the successful establishment of a rat model with SCI. After the completion of the experiment, the incision wound was washed with saline, sutured, and sterilized. Then, the rats were placed in an incubator until awakening. The surgical incisions of the rats with SCI were disinfected with iodine solution twice a day for three consecutive days. Meanwhile, each rat was subjected to intramuscular injection of penicillin sodium in the hind limb at a daily dose of 40,000 U for three days. The lower abdomen was massaged with a gentle squeezing of the bladder several times each day to help with urination and defecation.
The rats were divided into the sham operation group, the SCI model group, and the DG treatment group, with 40 rats in each group. Each group was divided into four subgroups according to one, seven, 14 and 21 days (d) post-injury, respectively. The rats from the sham operation group were subjected to the exposure of the spinal cord without weight-striking damage of the spinal cord. The DG group was subjected to gavage with DG solution at a dose of 100 mg/kg, while the SCI and sham groups were subjected to the gavage with normal saline at an identical volume. The gavage was conducted according to the required time points once a day for each group.

Functional Evaluation
After SCI, behavioral testing for motor function was performed once a day for 21 days. Functional recovery of each rat was evaluated according to the Basso, Beattie, and Bresnahan (BBB) scale [51], which is composed of 21 different criteria for the movement of the hindlimb (from no observable hindlimb movements to normal locomotion). This scale is based on precise observation of hindlimb movements, stepping, and coordination for 4 min in an open field by a trained observer who was blind to experimental conditions.

Sampling
After anesthesia, five rats in each group were subjected to aorta perfusion with 0.9% saline. After the effluent was clear, perfusion with 400 mL of 4% neutral paraformaldehyde for 60 min was conducted. Then, the T9-T11 spinal cord was harvested for fixation with through 4% neutral paraformaldehyde and stored for future use. The rest of the rats within the group were perfused with 0.9% saline through the aorta. After the effluent was clear, the T9-T11 spinal cord was removed, rapidly frozen in liquid nitrogen, and then stored in a refrigerator at −80 • C for future use.

Hematoxylin-Eosin (HE) and Nissl Staining
Fixed spinal cord tissue was subjected to dehydration, dipping in wax, embedding, and sectioning at the thickness of 5 µm in preparation for parallel HE and Nissl staining. Continuous coronal sectioning was performed, and one slice was taken every 10 slices. The change in structure and neuronal Nissl bodies in the spinal cord was observed under a common optical microscope.

Immunochemistry
The wax-embedded spinal tissue in each group was sectioned into slices at a thickness of 5 µm. After wax removal, hydration and antigen repairing, the section was rinsed with PBS (phosphate buffer saline) for three times for 3 min each time. Then, the section was subjected to blocking using 10% BSA (bovine serum albumin) for 30 min. After washing with PBS 3 times for 3 min each time, the primary antibody at a dilution of 1:250 was added and incubated overnight at 4 • C. After washing with PBS three times for 3 min each time, the secondary antibodies for FITC and CY3 at a dilution of 1:1000 were incubated with the tissue slices at room temperature for 1 hr. The slices were washed with PBS 3 times for 3 min each time. DAPI staining and anti-quenching and sealing were conducted. The stained slices were subjected to the observation under a fluorescence microscope. In total, five slices were collected from each group, and three different visual fields were selected for image acquisition under a fluorescence microscope at a magnification of 400×. The slice incubated with PBS instead of primary antibody was used as the negative control. The analysis of the mean optical density (MOD) for neurons was performed using Image-pro plus 6.0 software (plus 6.0, Media Cybernetics, Bethesda, MD, USA).

Transmission Electron Microscopy Evaluation
The rats from different groups were sacrificed at one, three, seven, 14, and 21 days after insult. The tissue samples were collected at a point of 5 mm T9-11 spinal ventral horn. Following fixation in 2.5% (w/v) glutaraldehyde overnight, spinal cord tissues were post-fixed in 2% (v/v) osmium tetroxide and blocked with 2% (v/v) uranyl acetate. Following dehydration in a series of acetone washes, tissues were embedded in Araldite for coronal sections. Semi-thin sections and toluidine blue staining were performed for observation of location. Finally, ultra-thin sections of at least three blocks per sample were cut and observed using a Hitachi-7500 TEM (Tokyo, Japan).

Terminal-Deoxynucleotidyl Transferase (TdT) Mediated Nick End Labeling (TUNEL) Assay
To enable TUNEL labeling to better illustrate neuronal apoptosis of neurons, TUNEL/DAPI double labeling was performed. TUNEL staining was performed using the In Situ Cell Death Detection Kit (Beyotime, Shanghai, China. C1088) according to the manufacturer's instructions. The number of positive cells was counted in 10 randomly chosen fields within each section. The sections per animal were respectively selected for incubation with DAPI for 3 min prior to the visualization of the sections with an optical microscope. The positive cells were identified by green TUNEL dots located in neurons with a blue nucleus at high magnification. The number of apoptotic neurons was calculated.

Protein Extraction from Spinal Cord Tissue and Western Blot
The spinal cord tissue was washed twice with PBS. RIPA (Radio Immunoprecipitation Assay) lysis buffer, protease inhibitor, and PMSF (Phenylmethanesulfonyl fluoride) were added to the sample. Then, the spinal cord tissue was subjected to homogenization on ice and centrifugation at 4 • C at 12,000 rpm for 15 min. Protein concentration was determined through the bicinchoninic acid (BCA) method. The proteins in the spinal cord tissue were denatured at 100 • C for 5 min. After cooling, the proteins in the samples were subjected to SDS-PAGE (Sodium dodecyl sulfate-polyacrylamide gel electrophoresis), transferred to a PVDF (Poly(vinylidene fluoride)) membrane at 275 mA for 90 min, blocked with 5% skim milk for 2 ht, and subjected to TBST (a mixture of tris-buffered saline (TBS) and Polysorbate 20 (also known as Tween 20)) washing three times for 5 min each time. Then, probing was performed with primary antibody of anti-Rheb (1:500), Bcl-2, Bax, caspase-3, AKT, p-AKT, mTOR, p-mTOR, LC3, p62, p70S6K, and β-actin (1:1000) with incubation overnight at 4 • C. Secondary antibody of HRP(Horseradish Peroxidase)-labeled goat anti-rabbit IgG(Immunoglobulin G) or HRP-labeled goat anti-mouse IgG was incubated with the membrane at room temperature for 4 hr. The images of the target proteins were developed, collected, and analyzed by the LAS4000 mini gel-imaging system.

Real-Time PCR of miR-155-3p
The total RNA and miRNA of collected T9-T11 spinal cord was extracted using a TRIzol extraction kit (Invitrogen, California, Calif., USA) and miRNA purification kit (Tiangen, Beijing, China) according to the manufacturer's instructions. Primers were designed using the software Primer 5 based on the sequence information from the miRBase database (http://www.mirbase.org) and GenBank (https://www.ncbi.nlm.nih.gov/genbank/), and synthesized by Biotechnology Engineering Co., Ltd. (Shanghai, China).
The miRNA-155-3p gene primer sequence was 5 -CCTCCTACCTGTTAGCATTAAC-3 . The U6 small nuclear RNA (snRNA) was selected as the internal reference for miRNA reverse transcription and qRT-PCR. The entire amplification reaction system was performed according to the instructions of the miRNA Real-Time PCR Assay Kit (SuperReal PreMix Plus (SYBR GREEN)) (Tiangen, Beijing, China). Reverse transcription of miRNAs was conducted using cDNA synthesis of miRNA first-strand (tail addition method). The Poly(A) tail addition reaction and the cDNA synthesis reaction were carried out simultaneously (reaction system of 20 µL). In brief, 2000 ng of total RNA, 10 µL of 2× miRNA RT solution mix, and 2 µL of miRNA RT Enzyme mix were added to RNase-free water to a final total volume of 20 µL. The reaction conditions were: bathing at 37 • C bath for 60 min, and heating at 85 • C for 5 min for enzyme inactivation, and preservation at 4 • C. The obtained cDNA reaction solution was diluted 20-fold to use as the template for quantitative fluorescence detection. The qPCR reaction system (20 µL of total reaction system) contained 10 µL of 2× SuperReal PreMix Plus (with SYBR), 0.6 µL of miRNA primer Universal PCR Primer R (miR155-3P) (10 µmol/L), 0.6 µL of Universal U6 Primer F (10 µmol/L), 4 µL of miRNA first-strand cDNA obtained through reverse transcription, and 4.8 µL of RNase-free water. Reaction conditions were as follows: pre-denaturation at 95 • C for 15 min, denaturation at 95 • C for 10 s, annealing and extension at 60 • C for 30 s, and 40 reaction cycles. The samples were treated by PCR in triplicate. The relative expression of miRNA-155-3p was calculated by the 2 −∆∆Ct method. Ct = Ct (miRNA) − Ct (U6), Ct = Ct experiment − Ct control. Experimental data were automatically collected and analyzed by the Mx3005P real-time fluorescence quantitative system (Agilent, Foster, CA, USA).

Statistical Analysis
All data were expressed as mean ± standard deviation (M ± SD). Statistical analysis was performed with GraphPad Prism 5.0 software using a one-way analysis of variance followed by a Bonferroni post-hoc test. A statistically significant difference was considered at p < 0.05.

Conclusions
DG can activate autophagy, down-regulate miR-155, and activate the Rheb/mTOR signal pathway, thereby alleviating spinal cord injury, promoting functional recovery of the spinal cord, and providing a theoretical reference for the treatment of spinal cord injury and the development and utilization of natural products for health promotion and neurodegenerative disease treatment. However, its exact mechanisms still need to be further explored.


extfieldsofstudy: Biology Medicine

text:
Characterisation of the NRF2 transcriptional network and its response to chemical insult in primary human hepatocytes: implications for prediction of drug-induced liver injury

The transcription factor NRF2, governed by its repressor KEAP1, protects cells against oxidative stress. There is interest in modelling the NRF2 response to improve the prediction of clinical toxicities such as drug-induced liver injury (DILI). However, very little is known about the makeup of the NRF2 transcriptional network and its response to chemical perturbation in primary human hepatocytes (PHH), which are often used as a translational model for investigating DILI. Here, microarray analysis identified 108 transcripts (including several putative novel NRF2-regulated genes) that were both downregulated by siRNA targeting NRF2 and upregulated by siRNA targeting KEAP1 in PHH. Applying weighted gene co-expression network analysis (WGCNA) to transcriptomic data from the Open TG-GATES toxicogenomics repository (representing PHH exposed to 158 compounds) revealed four co-expressed gene sets or ‘modules’ enriched for these and other NRF2-associated genes. By classifying the 158 TG-GATES compounds based on published evidence, and employing the four modules as network perturbation metrics, we found that the activation of NRF2 is a very good indicator of the intrinsic biochemical reactivity of a compound (i.e. its propensity to cause direct chemical stress), with relatively high sensitivity, specificity, accuracy and positive/negative predictive values. We also found that NRF2 activation has lower sensitivity for the prediction of clinical DILI risk, although relatively high specificity and positive predictive values indicate that false positive detection rates are likely to be low in this setting. Underpinned by our comprehensive analysis, activation of the NRF2 network is one of several mechanism-based components that can be incorporated into holistic systems toxicology models to improve mechanistic understanding and preclinical prediction of DILI in man. Electronic supplementary material The online version of this article (10.1007/s00204-018-2354-1) contains supplementary material, which is available to authorized users.

Introduction
Drug-induced liver injury (DILI) remains a leading cause of clinical liver failure and a hindrance to the development of new therapies, making it a priority area for the advancement of new risk prediction methods (Chen et al. 2014). Hepatocellular perturbation is often coupled to the activation of stress-responsive transcription factors that regulate networks of cytoprotective genes. By reflecting activation of stress-inducible transcriptional networks, such genes are considered useful mechanistic markers of the cellular consequences of chemical insult (Wink et al. 2014). The increasing availability of public data sets describing transcriptional responses to large numbers of compounds provides an excellent opportunity to explore the behaviour of key stress-responsive networks from a broad perspective. One example of such a data set is the Japanese Toxicogenomics Project's Open TG-GATES repository, which contains transcriptomic data from primary human hepatocytes (PHH) exposed to 158 compounds (primarily therapeutic drugs, but also experimental toxins, cytokines and growth factors) at up to three different concentrations and time points (Igarashi et al. 2015). This and other resources can support a deeper understanding of the molecular make-up of stress-inducible transcriptional networks, which is key to the design of, and interpretation of data from, preclinical assays and in silico models that can reliably reflect their perturbation, and thus provide a necessary level of confidence in the understanding and prediction of relevant human risk (Oshida et al. 2015a, b, c;Sutherland et al. 2018).
Chemical and oxidative stresses are associated with the adverse effects of many compounds, particularly electrophiles and free radicals that react with critical macromolecules and disrupt normal redox processes (Park et al. 2011;Pereira et al. 2012). In mammalian cells, these insults are counteracted by an antioxidant stress response directed by the transcription factor nuclear factor erythroid 2 related factor 2 (NRF2 or NFE2L2) (Bryan et al. 2013). Governed by its repressor Kelch-like ECH-associated protein 1 (KEAP1), NRF2 regulates the basal and inducible expression of cytoprotective genes which contain antioxidant response elements (AREs) in their promoter and/or enhancer regions (Bryan et al. 2013). Despite the mechanistic association of NRF2 with DILI and other forms of drug toxicity in rodents (Clarke et al. 2016), and the increasing interest in using translationally relevant human platforms to improve hazard prediction in man, relatively little is known about the makeup of the NRF2 transcriptional network and its response to chemical perturbation in PHH, which are often used as a translational in vitro model for investigating DILI. To address this knowledge gap, we have used siRNA to modulate NRF2 activity and identify genes that respond robustly to perturbation of the network in freshly-isolated PHH. Using an unbiased weighted gene co-expression network analysis (WGCNA) approach (Langfelder and Horvath 2007;Sutherland et al. 2016) we have also exploited the wealth of data available in TG-GATES to highlight the complex response of the NRF2 network to a range of chemical insults in PHH, and considered how perturbation of this pathway can be used in the context of understanding why certain compounds are associated with DILI in man. Our findings have important implications for the ongoing development of in vitro platforms and in silico models designed to improve the prediction of clinical DILI risk at an early stage of preclinical drug development.

PHH isolation, culture and drug treatment
Liver tissue was obtained from the Liver Cell Lab at the Karolinska University Hospital (Huddinge, Sweden) or Aintree University Hospital (Liverpool, UK) by qualified medical staff. All patients donated tissue as part of planned liver resections for various indications (Table S1). Informed consent in writing was obtained from each patient and the study protocol conformed to the ethical guidelines of the 1975 Declaration of Helsinki. Immediately following removal from the patient, excess healthy liver parenchyma was separated from the specimen and placed in cold Eagle's minimum essential medium, and transported to the laboratory on ice. Tissue dissociation and hepatocyte isolation were performed using a two-step collagenase perfusion procedure, as described previously (Strom et al. 1996). Only cell preparations with viability of ≥ 75% (as per Trypan blue exclusion test) were used for the experiments described here. The cell suspension was diluted in William's medium E without phenol red, supplemented with 25 mM HEPES and 2 mM l-glutamine, pH adjusted to 7.4 (modified William's medium E) plus 10% FBS. Cells were seeded at a density of 1 × 10 6 cells/mL onto Type I collagen-coated plates (Beckton Dickinson, UK) and maintained at 37 °C in a 5% CO 2 atmosphere. After 3 h, the medium was replaced with fresh modified William's medium E not supplemented with FBS. After being allowed to adhere to the plates for a total of 16 h, untransfected cells were exposed to 1 mM diethylmaleate for 2 h, or 10 µM sulforaphane for 24 h. Compounds were dissolved in DMSO and the solvent content of the media was controlled to 0.5% in all cases.

siRNA transfection
For siRNA transfections, siRNA duplexes targeted against either human NRF2 (D-003755-05, subsequently referred to as siNRF2) or human KEAP1 (D-012453-03, subsequently referred to as siKEAP1), and a scrambled non-targeting control siRNA duplex (D-001210-03, subsequently referred to as siCON), were obtained from the Dharmacon siGENOME library (Thermo Fisher Scientific, UK). Immediately prior to plating, cells from individual donors were reverse-transfected with 20 nM siRNA using Lipofectamine RNAiMAX (Life Technologies, UK) in accordance with the manufacturer's instructions. Plated cells were maintained at 37 °C in a 5% CO 2 atmosphere for 48 h to enable NRF2 or KEAP1 knockdown.

Microarray analysis and bioinformatics
Total RNA (50 ng, free from genomic DNA) was labelled and amplified using a Low-Input Amplification Kit (Agilent, USA). Amplified Cy3-labelled RNA (600 ng) was fragmented and loaded onto SurePrint G3 Human Gene Expression 8 × 60K v2 arrays (Agilent). Following overnight hybridisation at 65 °C, the arrays were analysed at the Liverpool Centre for Genomic Research, according to the manufacturer's instructions, using an Agilent G2505C Microarray Scanner. The data were extracted using Agilent Feature Extraction software v11.0.1.1. Differential gene expression analysis was conducted using the limma package within the R programming environment (R-Development-Core-Team 2005), enabling simultaneous comparisons between multiple treatments using design and contrast matrices via a linear regression model. To account for inter-individual differences in basal gene expression, the level of each gene was determined in siNRF2-or siKEAP1-transfected cells relative to siCON-transfected cells derived from the same donor. The significance (raw P-value) of estimated log 2 fold changes for the contrasts was evaluated using limma function eBayes, and the impact of multiple testing was adjusted using the Benjamini and Hochberg approach (Benjamini and Hochberg 1995). Differentially expressed genes were defined as those with an adjusted P value < 0.05 when compared with the level in siCON-transfected cells. Ingenuity Pathway Analysis (IPA; http://www.ingen uity.com) enrichment statistics were used to reveal biological pathways perturbed in siRNA-transfected cells. Pathways represented by a single gene/protein were excluded for robustness. Gene ontology (GO) term enrichment analysis was performed using GOrilla (cbl-gorilla.cs.technion.ac.il).

WGCNA and bioinformatics
Affymetrix HGU133-2 microarray CEL files generated from all PHH experiments were downloaded from the Open TG-GATES repository, jointly normalized using Robust Multiarray Average (RMA) using the Affy R package. Brainarray CDF (version 19) annotation were used to map probe sets to Entrez IDs (http://brain array .mbni.med.umich .edu/ Brain array /Datab ase/Custo mCDF/genom ic_curat ed_CDF. asp). Under this annotation, every gene is defined by a single probe set. This resulted in 17,500 probe sets, each mapping to a single gene, used for analysis. The TG-GATES repository contains 941 PHH experiments. An experiment denotes expression results for PHH treated with a given combination of compound, concentration and time compared to time-matched PHH treated with DMSO. For each experiment, log 2 fold change values were calculated for all genes by subtracting average log2 intensity for DMSO arrays from average log 2 intensity for treatment arrays. To identify co-expressed genes from the PHH data, we used the WGCNA R package (Zhang and Horvath 2005) and applied it to a matrix consisting of 941 rows (PHH experiments) and 17,500 columns (log 2 fold change values for probes). We created unsigned modules (i.e. grouping together co-induced and repressed genes). The soft power-selection algorithm in WGCNA produces asymptotic curves which requires selection of an arbitrary level of agreement with a scale-free network topology (e.g. 90%). Small changes in the threshold can lead to large changes in the selected power. We reasoned that non-expressed genes in PHH, which represent 'noise' in microarray experiments, should not be co-expressed and, therefore, not contained within any module. Therefore, we created modules for each soft-power parameter setting of 4, 6, 8 and 10, and performed t tests assuming unequal variance on the category "yes/no" (indicating whether a gene is member of any module) vs. genes' intensity in DMSO control. We selected 6 as the optimal soft-power parameter, as it maximized the t-statistic. Using WGCNA to merge similar modules (having correlation of their eigengene values ≥ 0.8), we obtained 399 modules containing 10,275 genes. As described previously (Sutherland et al. 2016(Sutherland et al. , 2018, for each experiment we calculated the eigengene (or module score) which summarizes log 2 fold change of their constituent genes. Briefly, this protocol deviates from the standard WGCNA approach in two steps: (a) when calculating a module score for a given experiment, log 2 fold change values are scaled only (not centred and scaled) to avoid producing non-zero module scores for experiments where all underlying genes are unperturbed, and (b) the raw module score is normalized to unit variance (Z-score scaling) facilitating comparison across modules and across treatments. A normalised module score summarises the level of gene activation or repression caused by a given treatment in the context of the large collection of compound perturbations. As such, an eigengene greater than + 2.0 or smaller than − 2.0 is a large (and significant) perturbation in the context of the 941 experiments. Finally, we calculated the correlation between a module eigengene vs. underlying genes' log 2 fold change across the 941 experiments, indicating the extent to which the eigengene summarizes the underlying gene's behaviour. To identify WGCNA modules that were enriched for genes found in the siRNA screen, canonical hypergeometric tests and subsequent correction for multiple comparisons (Bonferroni) were performed amongst modules that contained at least one gene identified in the siRNA screen.

Calculation of performance indicators
An eigengene value of ≥ 2.0 was considered a positive (and significant) perturbation of a module. When a compound provoked a positive response at several time points and/or concentrations, the experiment with the highest eigengene value was used. For the combined metric, an eigengene value of ≥ 2.0 for any of the NRF2-associated modules was considered a positive perturbation. Negative control modules were selected using the following criteria: (a) comparable number of genes (6-11) to the NRF2-associated modules, (b) absence of enriched GO biological processes associated with the response to oxidative stress and related terms. After classifying the compounds based on features associated with DILI, as described in the supplementary methods, performance indicators were calculated as follows (using compounds with intrinsic biochemical reactivity as an example): sensitivity, positive perturbations by intrinsically reactive compounds divided by the total number of such compounds; specificity, 1 minus (positive perturbations by non-reactive compounds divided by the total number of such compounds); accuracy, (positive perturbations by intrinsically reactive compounds plus non-perturbations by non-reactive compounds) divided by the total number of intrinsically reactive and non-reactive compounds; positive predictive value, positive perturbations by intrinsically reactive compounds divided by the total number of positive perturbations (i.e. by intrinsically reactive and non-reactive compounds); negative predictive value, non-perturbations by non-reactive compounds divided by the total number of non-perturbations (i.e. by intrinsically reactive and nonreactive compounds). ROC curve analysis was performed using GraphPad Prism version 7.

Supplementary methods
Details of the methods used for quantitative PCR, determination of cellular glutathione content, immunoblotting and the classification of compounds based on features associated with DILI are provided as supplementary material.

Genetic modulation of NRF2 activity in PHH
In light of our interest in characterising the NRF2 transcriptional network in PHH, we first examined our previously published proteomic data sets (Bell et al. 2016;Heslop et al. 2017) to confirm that canonical NRF2 targets are expressed at comparable levels in PHH and the liver tissue from which they were isolated, and are relatively stable in cells cultured in 2D for up to 7 days, in contrast to e.g. cytochrome P450 drug metabolising enzymes which are known to exhibit a marked time-dependent decrease in expression under culture conditions (Fig. S1). These analyses justified the use of cultured PHH for investigating NRF2-regulated processes. Next, we sought to modulate NRF2 activity in PHH from four donors (Table S1A) by transfecting them with either siRNA targeting the transcription factor (siNRF2), siRNA targeting KEAP1 (siKEAP1) or a scrambled, non-targeting control duplex (siCON). Knockdown of NRF2 and KEAP1 mRNA in cells from each donor was confirmed by qPCR ( Fig. 1a, b) and immunoblotting (Fig. 1d), and was associated with a decrease and increase, respectively, in expression of the canonical NRF2 target gene NQO1 at both mRNA and protein levels (Fig. 1c, d). In addition, significant changes in glutathione content were detected in PHH transfected with siNRF2 and siKEAP1, compared with the level in cells transfected with siCON ( Fig. 1e). These data demonstrate the successful genetic modulation of NRF2 activity in PHH and the expected impact on thiol-based redox status.

Characterisation of the NRF2 transcriptional network in PHH
To provide a detailed characterisation of the NRF2 network in PHH, we performed microarray analysis on RNA extracted from the cells transfected with siCON, siNRF2 or  (Table S3). The response of NQO1 in the microarray analysis of cells transfected with siNRF2 or siKEAP1 showed good correlation with the targeted qPCR analysis performed on the  (Table S1A) were transfected with 20 nM siCON, siNRF2 or siKEAP1 for 48 h. a-c qPCR analysis of a NRF2, b KEAP1 and c NQO1 mRNA levels in cells transfected with siCON, siNRF2 or siKEAP1. mRNA levels are normalised to GAPDH and expressed as a percentage of the levels in siCON-transfected cells. d Immunoblot analysis of NRF2, KEAP1 and NQO1 protein levels in cells transfected with siCON, siNRF2 or siKEAP1. Representa-tive blots from one PHH donor are shown. β-Actin was probed as a loading control. *Non-specific antibody signals. e Quantification of GSH levels in cells transfected with siCON, siNRF2 or siKEAP1, or exposed to 1 mM diethylmaleate (DEM) for 2 h. Data represent mean + SD. of PHH from n = 4 donors, transfected and analysed separately. Statistical analysis of qPCR and GSH data was performed with a paired t test; *P ≤ 0.05 , ***P ≤ 0.001 versus siCON Microarray analysis was performed on PHH, isolated from four donors, 48 h after transfection with siCON, siNRF2 or siKEAP1. a, b Volcano plots depicting differentially expressed genes in cells transfected with a siNRF2 and b siKEAP1, compared with cells from the same donor transfected with siCON. Each point represents a single gene probe, with those shaded orange (raw P value) and red (adjusted P value) found to be significantly different (P < 0.05) between groups. c Overview of microarray data sets and the strategy for the identification of genes that are highly sensitive to modulation of NRF2 activity. DE differentially expressed. The relevant supplementary tables are indicated in the boxes same samples (Fig. S2), indicating the robustness of the microarray data. IPA mapping of the 993 probes that were significantly down-regulated in siNRF2-transfected cells, or the 1176 probes that were significantly up-regulated in siKEAP1-transfected cells, demonstrated that the 'NRF2mediated Oxidative Stress Response' was one of the most significantly enriched pathways in both cases (Table S4), further supporting the genetic modulation of NRF2 activity in PHH. Other pathways commonly affected by both siNRF2 and siKEAP1 included 'Xenobiotic Metabolism Signaling', 'Glutathione Mediated Detoxification' and the 'Pentose Phosphate Pathway' (Table S4), consistent with established biological roles of NRF2. GO term enrichment analysis yielded comparable results (Table S5).

Fig. 3
Validation of NRF2-regulated genes in PHH. a Microarray, b qPCR and c immunoblot determination of the indicated genes/ proteins in PHH 48 h after transfection with siCON, siNRF2 or siKEAP1. a Note that the Agilent SurePrint G3 Human Gene Expression 8 × 60K v2 arrays contained a probe (A_23_P129903) that is homologous to sequences in the TRIM16 and TRIM16L genes. b By qPCR, TRIM16L-specific primers yielded data that essentially matched that generated in the microarray, whilst TRIM16 was barely detectable (data not shown). For qPCR, gene expression levels are normalised to GAPDH. c β-Actin was probed as a loading control. Representative blots from one PHH donor are shown. d qPCR determination of the indicated genes in a separate batch of PHH (Table S1B) 24 h after exposure to 10 µM sulforaphane. Data represent mean + SD of PHH from n = 4 donors, treated and analysed separately. Statistical analysis of qPCR data was performed with (a, b) a t test or d a Kruskal-Wallis (Conover-Inman) test; *P ≤ 0.05; **P ≤ 0.01; ***P ≤ 0.001 versus siCON or DMSO To identify genes that were most reflective of NRF2 modulation in PHH, and minimise the impact of potential off-target effects of the individual siRNAs, we focused on genes that were oppositely regulated in cells transfected with siNRF2 and siKEAP1. We identified 129 probes (representing 108 unique genes) that were significantly downregulated in siNRF2-transfected cells and upregulated in siKEAP1transfected cells and, therefore, judged to be positively controlled by NRF2 ( Fig. 2c; Table S6). Again, IPA mapping revealed that the 'NRF2-mediated Oxidative Stress Response' was amongst the most significantly altered pathways represented by these genes (Table S7), whilst 'oxidation-reduction process' was the most significantly enriched GO term (Table S8). In addition to genes known to be NRF2 targets in other cells (including AKR1B10, NQO1, PIR, SRXN1 and TXNRD1), our filtering strategy also revealed several putative novel NRF2-regulated genes, including F2RL2, LOC344887 and TRIM16L (Fig. 3a). We used qPCR and/or immunoblotting to confirm the altered expression of these genes in the siRNA-transfected PHH samples used to generate the microarray data (Fig. 3b, c). Furthermore, we demonstrated that these genes were upregulated in response to pharmacological activation of NRF2 by exposing PHH from four new donors (Table S1B) to sulforaphane for 24 h (Fig. 3d). Only 15 probes (representing 15 unique genes) were both significantly upregulated in siNRF2-transfected cells and downregulated in siKEAP1-transfected cells (Table S6), and thus judged to be negatively controlled by NRF2. IPA mapping identified 'Granulocyte Adhesion and Diapedesis' (related to the chemokines CXCL1 and CXCL2) as the only significantly altered pathway represented by the 15 genes (Table S7), consistent with GO term enrichment analysis (Table S8) and the reported link between NRF2 and the inhibition of nuclear factor kB and inflammatory processes (Wardyn et al. 2015). Taken together, these data highlight established and putative novel NRF2-regulated genes, representing important physiological pathways and functions, in PHH.

Mapping individual NRF2-regulated genes to networks using WGCNA
To further characterise the NRF2-driven gene network and provide a systems level insight into its response to chemical insult in PHH, we used WGCNA which provides a unsupervised and quantitative assessment of the interconnectedness of genes. Specifically, within a large set of gene expression data, WGCNA identifies genes that exhibit coregulated behaviour and organizes them into defined sets of highly interconnected genes or 'modules', thereby reducing dimensionality and providing an unbiased interpretation of network or pathway responses (Langfelder and Horvath 2007;Zhang and Horvath 2005). We first applied WGCNA to the entire TG-GATES PHH dataset, representing expression results for 17,500 probes in cells treated with one of 158 drugs/cytokines/growth factors at up to three concentrations and three time points, compared to time-matched PHH treated with DMSO (Grinberg et al. 2014;Igarashi et al. 2015). Using this approach, we identified 399 modules comprising 10,275 genes (full details to be published separately). Of the 108 NRF2-associated genes identified in our PHH siRNA screen (Fig. 2), 66 were found to be members of 31 modules (Table S9). Amongst the remaining 42 genes, 31 (including 17 unknown/uncharacterised transcripts) were not represented in the 17,500 probes within the TG-GATES data set. Therefore, 11 NRF2-associated genes that were present in the TG-GATES data set did not meet the co-expression criteria for module membership (Table S9). Whilst the 66 NRF2-associated genes were distributed across 31 modules, five modules (namely 2, 144, 192, 224 and 325) were significantly enriched for genes identified in the siRNA screen, suggesting tight co-regulation in response to NRF2 activation (Table 1). In light of the size of module 2 (743 genes, including 15 NRF2-associated genes identified in our siRNA screen) we focussed on the remaining four modules, which comprised 6-11 genes. The hub genes for these modules (those most correlated with the eigengene value, which summarises the induction or repression of the module as a whole (Langfelder and Horvath 2007)) were TXNRD1, CYP4F11, CD109 and CBR3, respectively. In addition, several of the remaining genes in the four modules were either down-regulated by siNRF2 or up-regulated by siKEAP1 in our PHH experiments (Table 1). Amongst these genes was F2RL2 (module 325), further highlighting it as a novel NRF2-regulated gene in PHH. In keeping with this, we noted that many of the genes in the four modules have been shown by ChIP-Seq to be direct NRF2 targets in human lymphoblastoid cells (Chorley et al. 2012), whilst the modules were represented by GO terms consistent with a role in xenobiotic metabolism and/or the oxidative stress response (Table 1 and S10). These observations support the enrichment of modules 144, 192, 224 and 325 with NRF2-associated genes that exhibit consistent co-regulated behaviour in PHH exposed to a broad range of chemical insults.

Response of NRF2-associated gene sets to chemical insult in PHH
Having identified modules that were enriched for NRF2associated genes with consistent co-expression behaviour across many treatments, we sought to distinguish compounds that produced significant perturbations in these modules, and by inference the NRF2 pathway, in PHH according to the TG-GATES data set. In total 47 compounds (including 34 therapeutic drugs) caused a significant change (i.e. eigengene value ≥ 2.0) of at least one of the four modules (Table 2). Notably, modules 144 and 192 responded to more compounds (36 and 35, respectively) than modules 224 (26 compounds) and 325 (19 compounds). Indeed, some compounds were found to cause the upregulation of selected NRF2-associated modules/genes but not others. For example, the established hepatotoxins azathioprine, diclofenac, flutamide and isoniazid provoked significant changes in four, three, two and one of the NRF2-associated modules, respectively (Fig. 4). These data highlight the complex responses of NRF2-associated modules and their constituent genes to chemical insult in PHH, and highlight the value of using modules (rather than selected individual genes) to quantify network perturbation.  (29) are indicated. Within each module, genes are ranked according to eigengene correlation, which indicates how well (hub-like) each individual gene is correlated with the module eigengene (1.0 is a perfect correlation). The top five enriched GO processes (ranked according to P-value) are shown for each module; the complete list for each module is provided in

Perturbation of NRF2-associated gene sets as an indicator of compound features associated with DILI
Finally, we sought to test whether the perturbation of NRF2associated modules, and by inference the stimulation of the NRF2 transcriptional network, could be a useful indicator of clinical DILI risk and/or compound features associated with hepatotoxicity. We surveyed the literature to classify the 158 TG-GATES compounds as; (a) those with intrinsic biochemical reactivity towards cellular biomolecules (10 compounds; Fig. 5a; Table S11), (b) those bioactivated to a chemically reactive metabolite in rodent and/or human in vitro or in vivo systems (83 compounds; Fig. 5b; Table S11), and (c) those are indicated. Bioactivation status and clinical DILI risk assignments as described in the main text  . 4 Response of NRF2-associated gene sets to chemical insult in PHH. Individual gene expression levels (% vehicle) and eigengene scores for modules 144,192,224 and 325 in PHH exposed to high concentrations of azathioprine (72.8 µM), diclofenac (400 µM), flutamide (50 µM) or isoniazid (10 mM) for the indicated times. All data were obtained from the TG-GATES PHH data set with experimental evidence for a lack of intrinsic biochemical reactivity and bioactivation (17 compounds; Table S11). To the best of our knowledge, the metabolic bioactivation status or intrinsic biochemical reactivity of the remaining 48 compounds (including seven cytokines and growth factors) has not been reported, and so these compounds were excluded from the relevant analyses. Separately, we found that 107 of the 158 compounds used in the TG-GATES PHH experiments were present in DILIrank (Chen et al. 2016), with 51 classified as 'most DILI concern' and 45 as 'less/ no DILI concern' (Fig. 5c; Table S12). The 11 compounds classified as 'ambiguous DILI concern', along with the 51 compounds not present in DILIrank, were excluded from the analyses.
In keeping with the responsiveness of NRF2 to direct chemical stress, perturbation of the NRF2-associated modules was found to be a very good indicator of intrinsic biochemical reactivity, with relatively high sensitivities, specificities and positive/negative predictive values in each case ( Fig. 5a and Fig. S3; Table S13). Importantly, the Fig. 5 Perturbation of NRF2-associated gene sets as an indicator of compound features associated with DILI. The 158 TG-GATES compounds were classified based on a intrinsic biochemical reactivity, b metabolic bioactivation status and c clinical DILI concern. Compounds classified as 'unknown' or 'ambiguous DILI concern' were excluded from the relevant analyses. Based on the module eigengene values for each compound (see Tables S12-15), performance indicators were calculated (see "Experimental Procedures" for details) to determine the association between module perturbation and the respective toxicity features. See Fig. S3 for performance indicators associated with individual NRF2-associated modules. P/NPV = positive/negative predictive value Table 3 Area under the curve (AUC) values for NRF2-associated gene sets as indicators of compound features associated with DILI Receiver operator characteristic (ROC) curve analysis was performed, based on the eigengene values of modules 144, 192, 224 and 325 for relevant compounds, to determine the association between module perturbation and the respective toxicity features (an AUC value of 1 = a perfect indicator, 0.5 = random). ROC curves are provided in Fig. S4 Module 144 Module 192  high sensitivity values associated with the perturbation of the NRF2-associated modules were not observed with a selection of similarly-sized gene sets (modules 140, 181 and 269; see Table S16 for compositions) enriched for GO processes not linked to the oxidative stress response ( Fig.  S4 and Table S13). NRF2 module perturbation was found to be a less sensitive indicator of metabolic bioactivation status and clinical DILI concern, with sensitivity values below 50% in all cases (Fig. 5b, c and Fig. S3; Table S14-15). However, specificity and positive predictive values remained relatively high in both contexts, suggesting that false positive detection rates should be low in these settings. These observations were supported by ROC curve analysis, with area under the curve (AUC) values above 0.8 for the prediction of intrinsic biochemical reactivity and generally below 0.6 for the prediction of metabolic bioactivation status or clinical DILI concern (Table 3 and Fig. S5). Using the two NRF2-associated modules (42 m and 320) identified in our previous WGCNA depiction of co-regulated gene sets in rat liver (Sutherland et al. 2016(Sutherland et al. , 2018, we observed similar trends (relatively low sensitivity but high specificity and positive predictive values for the prediction of metabolic bioactivation status or clinical DILI concern) when analysing the associated TG-GATES in vivo data set (Fig. S6), indicating that our findings are relevant across species and experimental models. Taken together, these analyses show that stimulation of the NRF2 transcriptional network is a very good indicator of intrinsic biochemical reactivity/chemical stress, but is a less sensitive indicator of metabolic bioactivation status and clinical DILI concern.

Discussion
There is an urgent need to develop a more holistic preclinical approach to relate chemical insults to relevant biological events to better predict DILI risk in man. Given its mechanistic association with hepatotoxicity in rodents (Clarke et al. 2016), NRF2 activation has been considered a relevant biological event that could be incorporated into next-generation models affording a greater understanding of DILI, and thus an improved prediction of clinical hazard. The selection of robust marker genes representing NRF2 activation, and other stress responsive transcriptional networks, is critical to such endeavours. In this study, we have used an unbiased approach, combining siRNA gene knockdown, transcriptomics and WGCNA, to reveal the molecular landscape of the NRF2 transcriptional network in PHH and identify co-expressed genes sets with consistent responses across many compounds. In doing so, we have determined the ability of NRF2 activation to reflect compound features associated with DILI. Our findings illustrate a number of important points. This study has demonstrated the high sensitivity of the NRF2-driven stress response in PHH to compounds with intrinsic biochemical reactivity towards cellular macromolecules, i.e. those compounds most likely to induce direct chemical stress. Whilst our findings will need to be confirmed with transcriptomics data from a larger selection of intrinsically reactive compounds, when such a resource becomes available, our observations are consistent with the canonical mechanism underpinning chemical stimulation of NRF2 signalling, whereby the transcription factor responds to agents that react with cysteine residues in KEAP1 to induce a chemopreventive phenotype at subtoxic concentrations (Bryan et al. 2013). Notably, all of the TG-GATES PHH experiments used minimally cytotoxic concentrations of each compound. Hence, the observed NRF2-driven responses can be regarded as transcriptional adaptations to relatively minor chemical insults, further exemplifying the sensitivity of the network to direct chemical stress.
It is important to note that chemical stress may or may not result in DILI and/or others forms of overt toxicity, depending on the extent of exposure and engagement of critical cellular targets. Consistent with this, we have shown that activation of NRF2 in PHH by itself is not a sensitive predictor of clinical DILI risk. However, in light of the complexity of DILI (reflected in the range of hepatotoxicity mechanisms represented by the TG-GATES compound set) it is inconceivable that a single cell signalling pathway could adequately reflect the range of chemical and biological perturbations necessary to afford a holistic prediction of DILI risk in man. Given that NRF2 has evolved as one of several stress-responsive transcription factors that serve to protect the liver and other organs from the consequences of a range of environmental burdens, it will be important to also consider the activities of these other cellular networks to better determine the balance between benign adaptation and deleterious progression in the face of a given chemical insult, thereby improving our mechanistic understanding and ability to predict clinical DILI risk. Such efforts will be aided by the emergence of large, granular data sets encompassing measurements of several stress response pathways alongside other pertinent endpoints (Huang et al. 2016;Wink et al. 2017Wink et al. , 2018. We have also shown the relatively low sensitivity of NRF2 network perturbation as an indicator of the propensity of a compound to undergo bioactivation to a reactive metabolite. Critically, reactive metabolite formation was not measured as part of the TG-GATES experiments, making it difficult to relate chemistry to biology in the same experimental setting. It is also known that reactive metabolite formation is neither necessary nor sufficient for all forms of DILI. Therefore, it is unlikely that all instances of drug bioactivation represent a genuine chemical insult, as exemplified by acetaminophen, which is safe at therapeutic doses despite 5-10% bioactivation to a reactive quinoneimine metabolite, and only toxic following overdose. NRF2 activation may, therefore, better reflect the cellular consequences of drug bioactivation, complementing chemistry-based reactive metabolite screens in the early phases of drug development, where consideration of drug exposure at relevant organs (e.g. through physiologically-based pharmacokinetic modelling) may support an improved understanding of the clinical relevance of stress response network perturbations observed in preclinical settings.
An interesting finding of our work is that the NRF2associated modules and their constituent genes respond in a complex manner to chemical insult. Indeed, we have identified several compounds that alter the expression level of some but not all NRF2 targets. Although the detailed underlying mechanisms require further investigation, the influence of other transcription factors, cellular signalling events and post-translational modifications may underlie the ability of a set of compounds to commonly activate the NRF2 stress response yet augment the expression of different subsets of NRF2 target genes, and/or to different degrees.
Using WGCNA, we have demonstrated that one approach to reducing this complexity is to consider genes as components of modules that, based on co-regulated behaviour across a large set of experiments, reflect perturbation of transcriptional networks associated with important cellular processes (Sutherland et al. 2016(Sutherland et al. , 2018. Such an approach can provide an unbiased means of selecting appropriate genes as the basis for targeted analysis of stress response network perturbations. Indeed, amongst the few published investigations of NRF2 signalling responses in PHH, most have used pre-curated lists (often based on other cell types) as the basis for the selection of NRF2-regulated genes (for example, see Souza et al. 2017). However, not all of these genes may be robust markers of NRF2 activation in PHH per se. In this regard, our demonstration that SRXN1 is oppositely regulated following siRNA knockdown of NRF2 and KEAP1 in PHH, and highly correlated with the eigengene of NRF2-associated module 144, supports our selection of this gene as a marker of NRF2 activation within a panel of HepG2 cell lines expressing green fluorescent proteintagged components of key cellular stress response pathways, which can be used to provide a deeper understanding of the consequences for cell health of activating NRF2 and other signalling pathways in the context of different chemical toxicities (Wink et al. 2017(Wink et al. , 2018. Interestingly, Srxn1 is not present in Nrf2-associated modules previously defined in rat liver (Sutherland et al. 2016(Sutherland et al. , 2018, further illustrating the utility of the module-based approach and context-specific nature of marker gene selection. We and others have previously used transcriptomics and proteomics to detail the biological processes that are altered in the livers of Nrf2-or Keap1-null transgenic mice (Kitteringham et al. 2010;Walsh et al. 2014;Wu et al. 2012). In keeping with the findings of those investigations, the present study has shown that, in PHH, NRF2 regulates the expression of a battery of genes that coordinate the response to chemical and oxidative stress, the disposition of xenobiotics, and the provision of NADPH and other cellular fuels, supporting a role for NRF2 in the maintenance of normal hepatic function in man. Our study also highlights several putative novel NRF2 target genes, including F2RL2, TRIM16L and the pseudogene LOC344887. Whilst further work is needed to define the biological significance of the interaction between NRF2 and these genes, LOC344887 (also known as NMRAL2P) was recently shown to respond to sulforaphane and regulate the NRF2-dependent induction of NQO1 in colon cancer cells (Johnson et al. 2017). In addition, our data highlight differences between the regulatory roles of NRF2 in the human and mouse liver, such as the high responsiveness to NRF2 modulation of members of the AKR family primarily in PHH, and the GST family predominantly in mouse liver (Kitteringham et al. 2010;Walsh et al. 2014;Wu et al. 2012). Therefore, whilst the over-arching cytoprotective function of NRF2 is conserved across mammals, it is clear that subtle differences in the regulatory roles of the transcription factor exist between species. It will be important to consider the preservation of the NRF2 network relationships described here when designing assays and in silico models intended to bridge across species. Examination of gene expression changes following noncanonical activation of the transcription factor via the glycogen synthase kinase-3 axis (Hayes et al. 2015) may further expand the landscape of the NRF2 network.
In summary, we have performed an unbiased characterisation of the NRF2 transcriptional network in PHH, used WGCNA to identify robust marker genes reflecting NRF2 activation, and shown that the perturbation of this signalling pathway is a very good indicator of the ability of a compound to cause direct chemical stress in cells. We have also shown that the activation of NRF2 signalling is a less sensitive indicator of clinical DILI concern, yet false positive detection rates are likely to be low in this setting due to the relatively high specificity and positive predictive values associated with perturbation of NRF2-associated gene sets by non-hepatotoxic compounds. In the early stages of preclinical drug toxicity assessment, high specificity assays are desirable such that safe compounds are not inappropriately terminated, whilst high sensitivity assays become more important at later stages as compounds move closer to human trials. Therefore, in this context, a positive NRF2 signal would likely merit further investigation of the toxic potential of a drug candidate. On the other hand, compounds that do not provoke a chemical stress response are unlikely to stimulate NRF2 signalling, although a lack of NRF2 signal could not guarantee safety. A better understanding of what NRF2 activation, and other relevant stress responses, can and cannot tell us about the risks associated with a given compound will contribute to the improved prediction of DILI, ensuring that drug development programs are not unnecessarily halted and supporting the progression of safe and effective new drugs into the clinic. Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creat iveco mmons .org/licen ses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.


s2fieldsofstudy: Biology Medicine

text:
Possible role of p53/Mieap‐regulated mitochondrial quality control as a tumor suppressor in human breast cancer

Mitochondria‐eating protein (Mieap), encoded by a p53‐target gene, plays an important role in mitochondrial quality control (MQC). Mieap has been reported to have a critical role in tumor suppression in colorectal cancer. Here, we investigated its role as a tumor suppressor in breast cancer. The enforced expression of exogenous Mieap in breast cancer cells induced caspase‐dependent apoptosis, with activation of both caspase‐3/7 and caspase‐9. Immunohistochemistry revealed endogenous Mieap in the cytoplasm in 24/75 (32%) invasive ductal carcinomas (IDC), 15/27 (55.6%) cases of ductal carcinoma in situ (DCIS) and 16/18 (88.9%) fibroadenomas (FA) (IDC vs DCIS; P = 0.0389, DCIS vs FA; P = 0.0234, IDC vs FA; P < 0.0001). In IDC, the Mieap promoter was methylated in 6/46 (13%) cases, whereas p53 was mutated in 6/46 (13%) cases. Therefore, the p53/Mieap‐regulated MQC pathway was inactivated in 12/46 IDC (26.1%). Interestingly, all tumors derived from the 12 patients with Mieap promoter methylation or p53 mutations pathologically exhibited more aggressive and malignant breast cancer phenotypes. Impairment of p53/Mieap‐regulated MQC pathway resulted in significantly shorter disease‐free survival (DFS) (P = 0.021), although p53 status is more prognostic in DFS than Mieap promoter methylation. These results indicate that p53/Mieap‐regulated MQC has a critical role in tumor suppression in breast cancer, possibly in part through mitochondrial apoptotic pathway.

driver of tumor development even in breast cancer. 9 p53 is clinically very important not only because of its high mutation rate but also because mutation is associated with more aggressive disease and worse overall survival. 10 p53 is a transcription factor that activates the expression of various downstream genes in response to DNA damage. 11 The central functions of this protein in tumor suppression are cell cycle arrest, apoptosis, DNA repair and anti-angiogenesis. [12][13][14][15][16] In particular, apoptosis is so important a function for p53-related tumor suppression that p53 activates target genes, including Bax, Noxa, Puma, Apaf-1 and p53AIP1, in response to DNA damage by radiation, UV and oxidative stress. 11,13,17 Although the mechanisms of apoptosis induced by DNA damage have been clarified, mitochondria are a pivotal organ for apoptosis, where these apoptosisrelated proteins localize and play an important role in mitochondria through caspase activation. 18 Recently, mitochondrial quality control has been revealed to be a novel function of p53. This function is regulated by a novel p53-inducible protein called mitochondriaeating protein (Mieap). 19 and NIX at the mitochondrial outer membrane may play a critical role in the translocation of lysosomal proteins from the cytoplasm to the mitochondrial matrix. 21 Although Mieap has been shown to be a key player in mitochondrial quality control, emerging evidence suggests that it also plays an important role in tumor suppression. In a mouse model, Mieap-deficient (adenomatous polyposis coli) ApcMin/+ mice had a much shorter lifespan and increased numbers and sizes of intestinal tumors compared to those in Mieap-wild type ApcMin/+ mice, suggesting that loss of Mieap increases ROS production in the intestinal mucosa and accelerates tumor progression. 22 Furthermore, Mieap-regulated mitochondrial quality control is frequently inactivated in human colorectal cancer by Mieap promoter methylation, BNIP3 promoter methylation or p53 mutation. 23 Considering that Mieap is a downstream target of p53, this novel mechanism for mitochondrial quality control is a new function of the p53 tumor suppressor. These results strongly prompted us to speculate that p53/Mieap-regulated mitochondrial quality control might be involved in tumor suppression in human breast cancer. In this study, we investigated this role.

| Fluorescence activated cell sorting analysis
We overexpressed Mieap using the constructed adenovirus in breast cancer cell lines. In addition, we treated these cells with etoposide (10 μmol/L). MCF-7, SK-BR-3 and MDA-MB-231 cells were infected with the adenovirus at a moi of 30 and 60, respectively; cells were then trypsinized and washed with PBS.
Subsequently, the cells were fixed in cold 70% ethanol at 4°C.
Fixed cells were washed and then incubated for 30 minutes at 37°C with PBS containing 0.2 mg/mL RNase, and then incubated for a further 30 minutes at approximately 22°C with PBS containing 20 μg/mL propidium iodide. Cells were analyzed using the EC800 Flow Cytometry Analyzer (Sony, Tokyo, Japan). For each analysis, 5000-20 000 cells were recorded, and the percentages of cells in different cell cycle phases (subG1, G1, S and G2/M) were determined using EC800 analysis software. The subG1 fraction was regarded as the apoptotic cell population. 25

| Caspase activity assay
Breast cancer cell lines were treated using the same conditions as described previously herein. Caspase activities were measured as described previously. 25 The cell pellets were resuspended in cell lysis buffer (50 mmol/L HEPES, 100 mmol/L NaCl, 0.1% CHAPS, 1 mmol/L EDTA, 10% glycerol, 0.1% NP40, pH 7.4) for 15-30 minutes on ice. Total protein (10-30 μg) was added to the reaction buffer (the same as the cell lysis buffer (except for the addition of 0.1% NP40) containing NAc-Asp-Glu-Val-Asp-pNA (Ac-DEVE-pNA, a substrate of caspase-3/7; Sigma-Aldrich), followed by incubation for 1-2 hours at 37°C. The free pNA cleaved from its precursor was quantified using a spectrometer at 405 nm (Bio-Rad, Hercules, CA, USA). The Caspase-Glo 9 Assay reagent (Promega, Madison, MI, USA) was used as the substrate for activated caspase-9. An equal volume of cell lysate and caspase-Glo 9 reagent was added to each well, mixed, and incubated in the dark for 1 hour at room temperature. The resultant luminescence, which is proportional to caspase-9 activity, was measured using a luminometer (Bio-Rad). Absorbance was determined using an Envision 2104 multilabel reader (Perkin-Elmer, Waltham, MA, USA). Data were analyzed using Student's t-test. A P-value <0.05 was considered statistically significant.

| Western blotting
Caspase plays a central role in apoptosis, and the cleavage of poly ADP-ribose polymerase (PARP) facilitates cellular disassembly.
Cleaved PARP is a marker of cells undergoing apoptosis. 26 Breast cancer cell lines were treated using the same conditions as those described previously herein. Cells were resuspended in radioimmunoprecipitation assay (RIPA) buffer (Sigma-Aldrich) with protease inhibitors for 15-30 minutes on ice. Protein (8 μg) was subjected to SDS-PAGE and then transferred to PVDF membranes. The membranes were blocked for 1 hour and then incubated with respective primary antibodies overnight at 4°C. The anti-Mieap antibody (rabbit, 1:1000) was prepared as described previously. 19 Anti-

| Immunohistochemistry
We examined the expression of Mieap in vivo using surgical speci- were incubated with EnVision reagents for 1 hour at room temperature. 22 IHC evaluation was performed by experienced pathologists.
When we started IHC for human cancer tissues in vivo, we repeated preliminary experiments for appropriate conditions using preimmune serum. For subsequent runs, we did not use this due to its shortage. However, as we found endothelial cells and smooth muscle of blood vessels to be good positive controls, we show Mieapexpressing endothelial cells in Figure S1.
Informed consent was obtained from all patients who agreed to provide a surgical specimen. This study was approved by the central ethics committee of Gifu University.

| Methylation-specific PCR and p53 mutation search
Methylation-specific PCR was performed as described previously. 19 The p53 mutation search was performed as described previously. 23 Genomic PCR was performed, followed by DNA sequence analysis.

| Statistical analysis
Immunohistochemistry for Mieap was analyzed using the χ 2 -test.
Disease-free survival (DFS) was estimated using the Kaplan-Meier method, and comparisons between groups were performed using a 2-sided log-rank test. Results were considered significant at P < 0.05.  Figure 1D).

| Mieap-induced cell death is apoptosis
As a preliminary experiment, we treated 3 cell lines with etoposide.
We found cleaved caspase-3/7/9 with activation of mitochondrial apoptotic factors particularly in SK-BR-3 and MDA-MB-231 ( Figure   S3). 17,27 To investigate the association between Mieap-induced cell death and caspases, we determined the expression levels and activation of caspase-3/7 and caspase-9 by western blotting and caspase activity assays. We found that SK-BR-3 was the most sensitive cell line to Mieap-induced cell death ( Figure 1D). Forty-eight hours after Ad-Mieap infection (30 moi), not only caspase-3/7 but also caspase-9 was significantly elevated in SK-BR-3 cells (caspase-

| Mieap expression is downregulated in tumors from breast cancer patients with an aggressive and malignant phenotype
We next performed IHC using surgical specimens obtained from 120  (Table 1).

| The p53/Mieap-regulated mitochondrial quality control pathway is genetically and epigenetically inactivated in aggressive and malignant breast cancer
To investigate inactivation of the p53/Mieap-regulated MQC pathway, we performed MSP for Mieap, NIX and BNIP3 promoters, as well as a p53 mutation search, because NIX and BNIP3 have been shown to be co-factors for MQC. 21 We investigated 46 IDC, which had been  Table 2). Of 12 cases, 6 with methylated Mieap promoter were luminal B and the other 6 with mutant p53 were triple negative as follows: 4, HER2-rich; 1, luminal B; 1, whose subtypes were clinically more aggressive and malignant phenotypes in breast cancer (Table 3). 30,31 Interestingly, the group with inactivation of the p53/Mieapregulated MQC pathway showed shorter DFS than that without these genetic alterations (P = 0.021, Figure 5). Four of six recurrent cases showed inactivation in the p53/Mieap-regulated MQC pathway. Of 46 cases, we separately showed DFS curves for individuals with or without Mieap promoter methylation or those with or without p53 alterations ( Figure S4). As shown in Figure S4A, there was no significant difference between those with Mieap promoter methylation (n = 6) and those without (n = 40). Six cases with a methylated Mieap promoter showed no expression of Mieap by IHC.
However, as shown in Figure S4B, there was a significant difference between those with p53 mutations (n = 6) and those without (n = 40;

| D ISCUSS I ON
According to a previous report, we demonstrated the physiological functions related to MQC, which were classified mainly into 2 mechanisms.
One is referred to as MALM, whereas the other is MIV, which occurs when salvage of the damaged mitochondria appears to be impossible for recovery and survival. 20   TA B L E 3 (Continued)

ACK N OWLED G M ENTS
We thank Shizu K for technical support and Enya K, Hirata M, Iwata A, Mori K and Takano K for administrative assistance.

D I S C LO S U R E S TAT E M E N T
KY has received grants and personal fees from Taiho