PGRMC1 effects on metabolism, genomic mutation and CpG methylation imply crucial roles in animal biology and disease

Progesterone receptor membrane component 1 (PGRMC1) is often elevated in cancers, and exists in alternative states of phosphorylation. A motif centered on PGRMC1 Y180 was evolutionarily acquired concurrently with the embryological gastrulation organizer that orchestrates vertebrate tissue differentiation. Here, we show that mutagenic manipulation of PGRMC1 phosphorylation alters cell metabolism, genomic stability, and CpG methylation. Each of several mutants elicited distinct patterns of genomic CpG methylation. Mutation of S57A/Y180/S181A led to increased net hypermethylation, reminiscent of embryonic stem cells. Pathways enrichment analysis suggested modulation of processes related to animal cell differentiation status and tissue identity, as well as cell cycle control and ATM/ATR DNA damage repair regulation. We detected different genomic mutation rates in culture. A companion manuscript shows that these cell states dramatically affect protein abundances, cell and mitochondrial morphology, and glycolytic metabolism. We propose that PGRMC1 phosphorylation status modulates cellular plasticity mechanisms relevant to early embryological tissue differentiation.

PGRMC1 phosphorylation site mutations cause pronounced metabolic changes, genomic mutation rates, and altered genomic CpG methylation, without affecting progesteronedependent doxycycline resistance. Hypermethylation of a S57A/Y180F/S181A mutant resembles pluripotent stem cells.

Background
Progesterone (P4) receptor membrane component 1 (PGRMC1) is a cytochrome b 5 (cytb5) protein and the archetypal a member of the membrane-associated P4 receptor (MAPR) family. It has a plethora of reported functions [1][2][3][4]. Human PGRMC1 contains predicted binding site motifs for Src homology 2 (SH2) and Src homology 3 (SH3) domain-containing proteins, with several other phosphorylation sites at S57, T178 and S181 being thought to additionally regulate these sites [5]. The SH3 target motif adjacent to S57, and the SH2 target motifs containing Y139 and T178/Y180/S181 are mutually juxtaposed on the folded protein surface, forming a potential proximity-stimulated tripartite signaling platform [5,6]. Thus, PGRMC1 represents a potential integration point and effector of many cell signaling pathways responsible for growth and proliferation [2,3,6,7].
The PI3K/Akt pathway has been shown to be affected by over expression of PGRMC1 [8], or 3x hemagglutinin (HA)-tagged PGRMC1 (PGRMC1-HA) phosphorylation site mutants in MCF-7 cells [9]. A PGRMC1-HA S57A/ S181A double mutant (DM) was associated with differential cell survival to oxidative stress. Furthermore, PGRMC1 was shown to be differentially phosphorylated between estrogen receptor-positive and negative breast cancers [9]. This study aims to further explore the biology of these PGRMC1 phosphorylation sites.
In a companion publication designed to investigate the biology of the DM mutation [10], we found that stably expressed DM PGRMC1-HA induced dramatic effects in MIA PaCa-2 (MP) pancreatic cancer cells relative to cells expressing wild-type PGRMC1-HA (WT), including elevated PI3K/Akt activity, altered metabolism, mitochondrial morphology and function, actin cytoskeleton proteins, and elevated migration rates relative to WT cells. To test the hypothesis that the DM effects required Y180, further mutation of S57A/Y180F/S181A generated triple mutant (TM) cells (which differ to DM cells by just the phosphate acceptor of Y180). Importantly, S57 and Y180/S181 (mutated in this study) probably occupy conformationally flexible regions in solution, immediately N-and C-terminally flanking the cytb5/MAPR domain. These are hence adjacent to each other in the folded protein structure [5,11], such that these mutations are unlikely to induce large conformational changes in the folded cytb5/MAPR domain. Indeed, these residues all appeared outside the MAPR domain during evolution of the animal lineage leading to humans [12].
TM cells exhibited reversal of some DM effects, including lowering of PI3K/Akt and migration activity, indicating that the DM phenotype was reliant upon Y180, and also led to dramatically impaired mouse xenograft tumor formation. Pathways enrichment analysis of a proteomics study suggested broad effects that were suggestive of plastic cellular reprogramming [10]. In this study we were interested in exploring which known or newly discovered PGRMC1 functions are related to these phosphorylation events. In this context it is informative to consider when different PGRMC1 functions may have been acquired during evolution because evolutionary conservation is a strong indicator of functional importance.
PGRMC1 interacts with and regulates ferrochelatase, the mitochondrial rate-limiting enzyme in heme synthesis [13]. The related protein PGRMC2 has recently been reported to chaperone heme from mitochondria to the nucleus in adipocytes, where heme regulates the transcription of genes that encode mitochondrial proteins. Lack of PGRMC2 led to altered mitochondrial morphology and function [29], superficially reminiscent of the mitochondrial effects of PGRMC1 phosphorylation mutants in MIA PaCa-2 cells [10].
PGRMC1 is induced by a variety of agents causing DNA damage, including the topoisomerase II inhibitor doxorubicin [30]. The PGRMC1 yeast Dap1 (damage-associated protein 1) homologs are associated with DNA damage response [17,31]. Heme-binding and cyP450 interactions are thought to be required for DNA damage response [17,30], which seems to be an ancient function.
PGRMC1 is strongly associated with sterol biology including association with the Insig/SCAP/SREBP complex that senses cholesterol levels and leads to mevalonate pathway induction [32] and production of the first sterol, lanosterol (reviewed by: [6,16]), which are apparently ancient functions. CyP450-interactions [1] conspicuously include PGRMC1 regulation of the most conserved eukaryotic cyP450 (lanosterol 14-alpha demethylase, CYP51A) to modify lanosterol from yeast to mammals which leads to cholesterol synthesis, providing the substrate for subsequent vertebrate steroidogenesis [6,17,18]. The PGRMC1mediated conferral of responsiveness to steroids like P4 presumably evolved in animals [33]. A complex between PGRMC1, the low density lipoprotein (LDL) receptor (LDLR) and the sigma-2 receptor transmembrane protein 97 (TMEM97) regulates the rapid internalization of LDLR via the membrane trafficking function of PGRMC1 [34], which is clearly an animal invention related to intercellular lipid transport, but related to ancient PGRMC1 steroid biology functions [16].
Anti-apoptotic activity of PGRMC1, especially to treatment with doxorubicin, has been reported by several studies and is linked to PGRMC1-dependent P4 responsiveness [35][36][37]. In mammals it involves PGRMC1 association with the CYP2D6 and CYP3A4 cyP450 enzymes responsible for doxorubicin hydroxylation and inactivation, and this activity required the PGRMC1 heme-chelating residue Y113 [11]. However, it remains unclear whether heme-chelation or tyrosine phosphorylation (and perhaps membrane trafficking) of Y113 is involved [5]. Therefore PGRMC1dependent P4-responsiveness may involve the membrane trafficking of actual P4 receptors to the cell surface, or PGRMC1 itself may mediate its own P4 response (reviewed by [2,30]). PGRMC1 undeniably confers P4 responsiveness to cells of the reproductive [3] and the central nervous systems [4]. A D120G heme binding-incapacitated mutant exhibited reduced P4-responsiveness [37], implicating hemebinding requirement for the P4 response. While steroid hormone signaling appeared with eumetazoans [33], apoptosis is employed widely by animals but it is widespread in unicellular and multicellular eukaryotes. Some underlying mechanisms were indeed inherited from protomitochondrial bacterial contributions to eukaryogenesis [38]. So it is difficult to judge whether anti-apoptosis is an ancient or an animal phenomenon.
To assay function, pharmacologists preferentially employ specific inhibitors. The AG-205 small molecule inhibitor of PGRMC1 is frequently used to study the role of PGRMC1 in various cellular processes. We have recently shown that AG-205 treatment results in the failure of many proteins associated with the actin cytoskeleton to immunoprecipitate with PGRMC1 [39]. Note that AG-205 is often cited as a PGRMC1-specific inhibitor. While it certainly affects PGRMC1-dependent function, AG-205 was modelled in silico against the heme-binding site of a plant MAPR protein [40]. As noted previously [2], there is no evidence that it is any more specific to PGRMC1 than to PGRMC2, or other MAPR family members Neudesin or Neuferricin, or indeed that it does not exert effects involving non MAPR proteins. Claims related to "PGRMC1specific inhibitor AG-205" should be treated circumspectly. Notwithstanding, we address AG-205 responsiveness of PGRMC1 mutations here, which are likely to reflect ancient MAPR biology.
In summary of the above, PGRMC1 has been reported to be involved in a plethora of different functions, potentially important to several central aspects of eukaryotic biology, and which seem to have evolved at different times. At least some of these seem to be regulated by phosphorylation events that appeared in eumetazoan animals [12], whereas other functions seem to be ancient in eukaryotes. Almost nothing is known about many functions. It remains unclear which functions are inhibited by the inhibitor AG-205, and how these are related to the PGRMC1-specific P4 response. There is an urgent need to systematically stratify and separate PGRMC1 functions. Here, we sought to further characterize the cells described in the companion paper [10] in an effort to better understand the biology. We report no detectable influence on P4-or AG-205 responses, but we observed pronounced effects on metabolism, genomic stability and epigenetic plasticity, implicating the PGRMC1 Y180 module with critical newly identified eumetazoan roles.

PGRMC1 phosphorylation state-dependent metabolic differences
In order to assess metabolic differences, cells were subjected to preliminary hyperspectral autofluorescence imaging: a non-invasive analytical technique [41]. Distinct spectral signatures of endogenous auto-fluorescent molecules [42,43] can provide valuable information about intracellular metabolic state [44]. In a pilot study we showed that MP and DM cells exhibited detectably different metabolites [43]. As a first step towards demonstrating altered PGRMC1-dependent metabolic processes, eight informative spectral features that discriminated the cells were chosen ( Fig. 1a and Table S1). The cells in each condition could be clearly distinguished using the autofluorescence of their endogenous metabolites. A pairwise linear discrimination analysis using the same spectral features indicated that all cell conditions were highly significantly separated from each other on the basis of endogenous autofluorescence (Fig. 1a). Figure 1b shows autofluorescence differences in spectral channel 3, thought to reflect flavin emission [44,45]. All cells over-expressing PGRMC1-HA proteins exhibited reduced levels, with DM and TM the lowest (Fig. 1b). We also noted significant differences observed between cell lines for autofluorescent emission at 700 nm (Fig. S1A), in a wavelength range at which heme-containing proteins and other porphyrins contribute strongly to autofluorescence [46]. Note that PGRMC1 not only binds heme but is involved in heme synthesis [13]. Heme itself does not emit fluorescence but leads to fluorescence quenching, however, various heme-containing proteins are fluorescent [46]. Although the actual identities of these and fluorescent species contributing to many hyperspectral channels remain unknown, several unidentified parameters also significantly discriminated between the metabolites present in the different cell lines, such as e.g. the ratio of channel 3 to channel 12 ( Fig. S1B,C, Table S1). These results indicated that PGRMC1 phosphorylation status dramatically alters cell metabolic state. In particular, the difference between DM and TM cells is due to the oxygen acceptor atom of Y180.

PGRMC1 phosphorylation mutants did not affect P4 or AG-205 responses
We were interested in whether DNA mutation rates may have been related to the P4-dependent protection against cell death, or to the mechanism of AG-205-induced cell death. Consistent with previously reported PGRMC1dependent anti-mitotic effects of P4 in Ishikawa endometrial cancer cells [36] and MDA-MB-231 breast cancer cells [47], incubation of cells in 1 μM P4 retarded cell proliferation of all cells over-expressing a PGRMC1-HA protein (WT, DM or TM), but not of MP cells relative to non-P4-treated control cells (Fig. 2a). When cells pretreated with or without P4 for 1 h were co-incubated in the presence doxorubicin (Dox), in the absence of P4 then WT, DM and TM cells were more susceptible to Dox-induced death ( Fig. 2b-c, white data points and boxes), however these cells exhibited greater P4-dependent protection against Dox-induced death ( Fig. 2b-c, shaded data points and boxes). As reported previously for Ishikawa cells [36], the 1 h pre-incubation with P4 was required to observe these protective effects of P4 (data not shown). There was no significant protective effect of P4  on cell survival for MP cells (Fig. 2c). Notably, there were no significant differences in P4-dependent protection between WT, MP or TM cells. There was also no difference in the sensitivity of any of these cells to AG-205-induced cell death, with essentially superimposable dose response curves, and a half maximal lethal effective concentration (EC 50 ) of approximately 32 μM AG-205 (Fig. S2). Therefore, although we observe PGRMC1-HA-dependent response to P4, this response is apparently unaffected by the influence of our PGRMC1 phosphorylation mutants. We conclude that the dramatic effects of altered PGRMC1 phosphorylation status observed in this and the accompanying manuscript [10] are due to different and newly described functions of PGRMC1, unrelated to the mechanism of P4induced vitality or AG-205-induced death. However, our mutant effects may require processes that emanate from the heme-binding MAPR domain, since the wild-type MAPR domain sequence was present in all mutants.

PGRMC1 phosphorylation status affects cytoplasmic redox status
Changes in metabolism could be associated with altered states of oxidative environment. Naphthalimide-flavin redox sensor 1 (NpFR1) is a fluorophore analogous to mitochondrial-specific NpFR2 [10,48], except that it is localized to the cytoplasm [49]. To further explore metabolic state, cells were treated with NpFR1 and fluorescent state was assayed by flow cytometry to reveal two cell populations: one with less oxidized NpFR1 and one with more oxidized NpFR1 (Fig. 3a). The type of experiment shown in Fig. 3 has been performed on numerous occasions. The actual ratio of low to high abundance distribution for a cell type tends to vary between experiments on different days, however the inter-cell comparisons of cells in the low fluorescent population between cell types examined on a single day remains relatively constant (MP > DM > WT > TM). The bimodal fluorescent peaks of Fig. 3a represent cell populations with different levels of cytoplasmic oxidation/reduction which have not been further characterized. All PGRMC1-HA-expressing cells exhibited the majority of cells in a more oxidized peak. TM (9.4%) and WT (17.7%) cells both exhibited less cells in the lower oxidized fraction (Fig.  3a,b, black arrows). For TM cells, which had most cells in the more oxidizing cytoplasmic population, the degree of NpFR1 oxidation in that fraction was significantly lower than for all other cells (Fig. 3a,c, white arrow). Higher cytoplasmic oxidative levels in TM cells suggested potential for effects on genomic mutation rates, associated with Y180 function.

PGRMC1 phosphorylation status affects pathways associated with replication, cell cycle, and mitogenic signaling
In the accompanying publication [10] we presented only highly significantly differential WebGestalt proteomics pathways with Benjamini-Hochberg adjusted p (adjP) values < 0.001. Several pathways related to cell cycle control and DNA replication were also detected, albeit with weaker significance, including five proteins predicted to be associated with decreased levels of ATM/ATR DNA damage sensing and control in TM cells (Fig. 4a), which was interesting to us because of the DNA damage associations of PGRMC1 (see Background). Although we could not detect differences in cultured proliferation of any cell lines [10], reverse phase protein array (RPPA) analysis revealed elevated levels of activated MKK4, MEK and ERK, as well as phosphorylated retinoblastoma protein (Rb) in TM cells (Fig. 4b). This suggests that PGRMC1 Y180 may be somehow involved in regulating the G1 checkpoint, which may be related to the reported ability of PGRMC1 and/or PGRMC1 to lower the propensity of cells to enter the G1 and the cell cycle, rather than enter G0, in spontaneously immortalized granulosa cells [20]. However, this hypothesis requires further investigation.

PGRMC1 phosphorylation status affects genomic integrity
Altered cytoplasmic oxidative environment, potential perturbations of the ATM/ATR pathways, and possible G1 checkpoint features without detectable differences in proliferation rates [10] led us to explore genomic stability. We cultured one cell line each of MP, WT, DM and TM for 30 passages, and sequenced the genomes of each. The PGRMC1-HA constructs were expressed at similar levels after 30 passages (Fig. 5a). Relative to the reference human genome sequence, the order of average mutation rates (difference to reference) per chromosome were WT < DM < TM (not shown), which was also reflected in unique mutations per MB, as rates per chromosome (Fig. 5b). This indicates that manipulation of the phosphorylation status of PGRMC1 Y180 can affect genomic stability, which is of potentially considerable relevance to cancer biology.

PGRMC1 phosphorylation status regulates the NNMT/1-MNA pathway
To further characterize metabolic differences between cells ( Fig. 1) we performed a pilot metabolomics study, where we observed elevated levels of 1-methylnicotinamide (1-MNA) in some but not all experiments (not shown), which we nevertheless chose to explore. 1-MNA is produced by the enzyme nicotinamide-N-methyl transferase (NNMT), which we found to be elevated in TM cells (Fig. 6a). NNMT utilizes S-adenosylmethionine (SAM) as a methyl group donor to convert nicotinamide to 1-MNA [50] (Fig. 6b).
Liquid chromatography/mass spectrometry (LC/MS) quantification revealed lower levels of SAM in DM and TM cells (Fig. S4A). Anti-NNMT shRNA attenuated NNMT levels in TM cells more than two-fold (Fig. S4B), which led to expected reduced levels of 1-MNA (Fig. S4C).
Contrary to the hypothesis that increased NNMT activity would reduce SAM levels as proposed for pluripotent human embryonic stem cells (hESCs) [52], SAM levels were not significantly altered by NNMT attenuation (Fig. S4D). We also treated DM cells with 1-MNA, which caused a morphological conversion from predominantly rounded to predominantly elongated form (Fig. S4E), suggesting that NNMT activity may contribute to the morphological differences between elongated MP and WT cells on the one hand and rounded DM and TM cells on the other. Other cell-types showed no significant change in cell morphology upon 1-MNA treatment (not shown).

PGRMC1 phosphorylation status regulates genomic methylation
Because of the reported effects of NNMT activity on global methylation levels [52], we assayed genomic CpG methylation levels. Hierarchical clustering of samples displayed striking separation of the cell groups, suggestive of highly divergent patterns of methylation among the cell types. All independently derived stable cell lines of each PGRMC1 state clustered closely to each other, and distant from other cell conditions (Fig. 7a).
We observed a highly significant number of differentially methylated probes (Fig. 7b). In the TM/DM comparison (differ by Y180 state) we observed a shift towards hypermethylation that was not as pronounced in WT/MP (differ by transfection of pcDNA3-1_PGRMC1-HA and hygromycin stable selection) or DM/WT (differ by PGRMC1-HA S57/ S181 status). We then considered the percentage of differentially methylated probes which were hypermethylated in each comparison and noticed a clear increase in the progression from WT/MP, DM/WT to TM/DM (Fig. 7c).
While the WT/DM effect could be due to both PGRMC1-HA expression and/or the stable antibiotic selection process, the differences between DM/WT and TM/DM were due solely to PGRMC1 mutation status. The overall level of CpG methylation was higher in TM cells across all chromosomes (Fig. 7c). Together, the results suggest a global increase in net methylation in this system, being highest in TM cells. Pathways enrichment suggests altered cell cycle control proteins. a Abundances of proteins detected as significantly differentially abundant in indicated Pathways Commons pathways by enrichment analysis (Accompanying paper [10]) at the adjP < 5% significance level. (b, c, e-g) Average reverse phase protein array (RPPA) normalized fluorescent intensity (NFI) from the indicated antibodies (see methods) is plotted from 6 replicate measurements. NFI is normalized to protein content. Significance levels are * p < 0.05, ** p < 0.01, *** p < 0.001. Conventions follow those of the accompanying paper [10]. d The ratio of average NFI of the anti-pRB antibody to that of the RB antibody Significantly differentially methylated probes were next assessed by genomic context region based on the EPIC manifest designations for relation to CpG island (Island; Open Sea, North Shelf, North Shore, South Shore, and South Shelf) because DNA methylation levels are known to depend upon genomic context [54]. The percentage of significantly hypermethylated probes was similar for all PGRMC1-HA-mutation states (WT, DM, TM), which differed from MP (Fig. S7A, left). By contrast, there were marked differences in the distribution of hypomethylated probes due to PGRMC1-HA mutation status, with DM and especially TM cells exhibiting elevated levels of hypomethylation at CpG islands and both North and South Shore sites, and reduced levels of hypomethylation (increased methylation) in Open Sea regions (Fig. S7A, right). This is consistent with increased levels of Open Sea methylation (in excess to Island and Shore sites) leading to the net increased methylation of TM cells, with a minority of local demethylation sites in the vicinity of CpG Islands at activated genes. See also the volcano plot of Fig. 7b, TM/WT comparison, where positive fold change (logFC > 0) is in excess.
The increased level of hypomethylation at Islands and Shores in TM cells was not detectable between 5′ untranslated and 3′ untranslated regions of annotated coding genes (Fig. S5). A proportion of cell-type and promoter-associated enhancers exhibited reduced hypermethylation or increased hypomethylation in the TM/ DM comparison (Fig. S6). This reduced level of enhancer methylation, normally associated with gene activation, occurred in the background of bulk overall increased CpG methylation in TM cells. Since TM and DM cells differ only in the presence of Y180 of the exogenous PGRMC1 protein, we conclude that Y180 phosphorylation status appears to specifically influence the genomic methylation status of cell-type specific enhancers.
We performed Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways and Gene Ontology (GO) enrichment analysis of the methylomics results, using probes from Island or Shore genomics regions only. For KEGG we separated the data into separate hyper-and hypomethylated gene sets. The KEGG results are presented in Fig. S7B, and full GO and KEGG results are provided as File S1. Each mutation state induced specific unique suites of detected pathways. Most significantly differential KEGG pathways were detected by all three cell comparisons.
Across all comparisons, the most significant enriched GO terms included some associated with embryology and developmental processes (Table S2, File S1). KEGG results reflected the activities of signal pathways, cell-extracellular interactions, and cancer biology (Table S2, File S1). E.g., for the TM/DM cell comparison, which differs in only the presence of phosphate acceptor oxygen of PGRMC1-HA Y180, KEGG results indicated that the most hypomethylated pathway in TM cells was PI3K-Akt (path:hsa04151). The PI3/Akt pathway was also among the top 5 in the hypermethylated data set for that comparison, and was significantly differential in all three comparisons (File S1). This means that genes associated with PI3K/Akt activity exhibited altered methylation states. We observed no significantly different alteration of CpG methylation sites associated with the Akt gene itself (not shown). In the accompanying publication we predicted reduced PI3K/Akt activity in TM cells on the basis of proteomics pathways enrichment, where the pathway was also significantly differential for all cell comparisons, and we demonstrated reduced phosphorylation of Akt substrates in TM cells [10].
By way of illustrative example, we further explored the biology underyling the pathways analysis by considering the unique KEGG pathways identified for both hyper- and hypo-methylated genes from Fig. S7A. 13 hypermethylated and 7 hypo-methylated pathways were unique to the TM/DM comparison (Fig. S7B). Their identities are given in Table S3. Hyper-methylated pathways included processes related to proteolysis and bacterial infection as well as DNA base excision repair and mismatch repair involved in averting mutations. The latter is striking because we also detected protein abundances involved in pathways associated with DNA damage (p53 effectors, and ATM/ATR) to be reduced in TM cells by proteomics (Fig. 4), with accompanying elevated mutation rate in TM cells (Fig. 5), which is fully consistent with the methylomics pathways biology. These circumstances suggest strongly that the study reflects true biology related to PGRMC1 Y180 function. Whereas NNMT was more abundant in TM cells (Fig.  6a), CpG sites in its gene were hypermethylated (not shown) as were the sites in many coding genes (Fig. S5). Histone methylation status, which we have not assayed, could be influential in determining gene expression levels, or a hypomethylated enhancer may activate transcription despite CpG methylation of the immediate gene locus. Approximately 80% of the genes encoding the 243 differential proteomics heat map of the accompanying paper [10] were significantly differentially methylated between one of the comparisons (not shown), suggesting that epigenetics may contribute greatly to the phenotypes observed for these cell types.
Our most important finding is that all tested PGRMC1 conditions induced specific fingerprint patterns of genomic CpG methylation, identifying PGRMC1 as a major epigenetic regulator. In particular, PGRMC1 Y180 affected pathways predicted to be associated with embryology, seemingly by targeting the genomic methylation status of CpG sites in a subset of cell-type and promoter-associated enhancers. Vertebrate embryology relies heavily on epigenetic histone and genomic methylation to direct cell plasticity between generations and in response to stimuli to determine the differentiation status of cells [55]. Y180 (and Y139) appeared in animal evolution in the common ancestor of eumetazoans (cnidarians, including corals, jellyfish, etc., and bilaterally symmetrical animals), concomitantly with the gastrulation organizer that drives the embryological tissue differentiation processes [12]. PGRMC1 Y180 was therefore evolutionarily acquired before the appearance of chordate organs, and the embryological differentiation processes that give rise to them, which immediately suggests how PGRMC1 may affect cancer cell plasticity. Differential PGRMC1 phosphorylation status may therefore be associated with the altered prognosis of triple negative breast cancers [9]. The increase of hypermethylated CpG sites observed in TM cells (Fig. 7, Fig. S5) is indeed reminiscent of the characteristic increase in hypermethylation associated with both hESCs and induced pluripotent stem cells (iPSCs). The genome becomes increasingly hypomethylated during embryogenesis and somatic tissue differentiation [56]. It is as if the maintenance of vertebrate differentiation processes depended upon the acquisition of Y180 and the organizer by the eumetazoan common ancestor [12], and, upon disturbance with Y180 functions in TM cells, the genomic CpG methylation pattern reverts to a phenotype reminiscent of the undifferentiated state of that ancestral evolutionary stage: prior to the appearance of deuterostome and chordate tissues and body plan. This phenomenon clearly merits future investigation.
NNMT regulates pluripotency in naïve hESCs, where it was proposed to deplete cellular SAM levels and pleiotropically lead to reduced histone methylation, including some genes in the Wnt pathway whose concomitant elevated expression affected pluripotency [52]. PGRMC1 involvement in hESC pluripotency [57] suggests that it could be mechanistically involved in the epigenetic modulation that regulates hESC pluripotency [52]. Perhaps the most striking result we observed was that each PGRMC1 phosphorylation mutant condition elicited a distinctly unique pattern of genomic methylation (Fig. 7), consistent with PGRMC1 playing a central role(s) in regulating genomic methylation, with each phosphorylated isoform capable of dynamically engendering different end-points. Clustering was applied using 520,031 probes with a coefficient of variation > 0.10. Height is calculated using the UPGMA method in the Lumi package and equates to the average distance between two members of two groups. b Volcano plots of the WT/MP, DM/WT and TM/DM comparisons. The log2 fold change (FC) (x axis) is plotted relative to the log10 p-value (y axis) to visualize the magnitude of the most significant features. Probes with an adjusted p-value < 0.05 and absolute log FC > 1 are plotted in green. Probes with an absolute logFC > 1 that were not significant are plotted in orange. Probes with a logFC < 1 that were significant at adjusted p-value < 0.05 are plotted in blue. Black probes were not significant. c Significant probes in each differential methylation comparison were extracted and the fraction of hypermethylated or hypomethylated probes per chromosome calculated. Methylation levels are higher in TM cells across all chromosomes As a target of Wnt signaling, TCF/LEF appeared very early in animal evolution, being absent from the single-celled choanoflagellate sister group to metazoans, but present in sponges [58]. By the time the common ancestor of eumetazoans evolved, embryological gastrulation initiated from the animal pole under control of Wnt/β-catenin signaling, inducing an axial organizer whose vertebrate descendant is known as the Spemann-Mangold organizer in amphibians. This initiates the orchestrated series of events where specific gene programs are activated by transcription factors that trigger differentiation into specialized cell types [58][59][60][61].
These processes must represent evolutionary novel metazoan functions. The correlation of appearance of the Y139/Y180 combination and organizer function [12] is consistent with the involvement of PGRMC1 in this process. This new hypothesized PGRMC1 functionality would seem to be superimposed upon more ancient PGRMC1 functions, such as sterol synthesis and steroid responsiveness. In support of this hypothesis, despite dramatic differences in energy metabolism, cell morphology, and epigenetics, all of our mutants responded similarly to P4-dependent resistance to doxorubicininduced death, which requires PGRMC1::cyP450 interactions [11], as well as to AG-205 toxicity (Fig. 2).
PGRMC1 is required to maintain the pluripotency of embryonic stem cells by influencing the p53 and Wnt/β-catenin pathways [57]. In the accompanying publication [10] we show that TM cells exhibit lower levels of PI3K/Akt activity, and lowered levels of GSK3-β phosphorylation on an Akt target site. The canonical Wnt pathway is activated when β-catenin translocates to the nucleus to act as a transcription factor after activation by Wnt ligands. Active GSK3-β phosphorylates β-catenin, targeting it for polyubiquitination and proteasomal degradation [62], and Akt phosphorylation of GSK3-β inhibits GSK3-β activity [62].
While our cells are not embryonic, lack of Y180 phosphorylation and low PI3K/Akt activity could correlate with increased GSK3-β activity, and lower β-catenin levels which would attenuate Wnt signaling via TCF/LEF. PGRMC1 has been associated with negative regulation of TCF/LEF in hGL5 granulosa/luteal cells [63] and hESCs [57]. Future experiments are required to explore this concept further in this and other cell models, which should be prioritized.
The propensity for animals and their cells to age was accompanied in evolution by the loss of pluripotent stem cells in the adult: an evolutionary necessity to establish larger collections of clonally communal but functionally specialized cells with diverse differentiation states [64]. That entire biology was associated with the evolution of the tissue and cell-type heterogeneity of bilateral animals, which started with the evolutionary appearance of the organizer [64]. The possible involvement of PGRMC1/PGRMC2 phosphorylation in cell cycle control mechanisms governing entry to G 0 and senescence [20,21,24,65] may accordingly be related to an overarching function associated with the reduction of the immortality potential of primitive eukaryotes, and the necessary evolution of a systematics of differentiation status and replicative control in bilateral animals [64]. This hypothesis merits future examination.
In mining our proteomics data from the accompanying manuscript [10] for enzymes related to the NAD/ NNMT/1-MNA pathway (Fig. 6b), we noted that nicotinamide phosphoribosyltransferase (NAMPT) was one of the relatively few proteins less abundant in TM compared to both DM and WT (Fig. S3). NAMPT catalyzes a competing reaction to NNMT, whereby nicotinamide is converted to nicotinamide mononucleotide, which is then transformed into NAD + by nicotinamide/nicotinic acid mononucleotide adenyltransferase. NAMPT is the rate-limiting enzyme in this pathway, and regulates sirtuin activity by modulating NAD + levels [76]. This entire pathway is related to aging processes [77], which were acquired by bilaterian animals along with tissue differentiation and cell cycle control mechanisms to evade cancer [64]. Indeed, these are the very types of biology with which PGRMC1 function is here implicated: suggesting that PGRMC1 phosphorylation could be involved in the mechanism(s) of animal aging. Considering the long coevolutionary history of bilaterians and PGRMC1 tyrosine phosphorylation sites, this merits future study.
Our results imply that the ability to interchange between phosphorylated states at T178, Y180, S181 and S190 [5] will be critical in supplying correct PGRMC1 function at appropriate times for the cell cycle and state of cell differentiation during early embryology. Clearly, Y180 is important for PI3K/Akt activity [10], a major pathway associated with aging [78].
As a final technical note, S181 has been discussed in the literature as a consensus CK2 site [2]. It must be noted in this context that PGRMC1 phosphorylation at S181 was unaffected by knockout of CK2 kinase activity in C2C12 mouse myoblast cells [79], indicating that at least in those cells CK2 was not responsible for S181 phosphorylation. Willibald et al. [80] recently used an antibody against PGRMC1 phospho-S181 (pS181) to argue that S181 phosphorylation is essential in PGRMC1 activation, as well as discussing the presence of phosphorylated PGRMC1 in immunohistochemistry of tissue biopsies [81]. A caveat to these results is that the peptide epitope used to generate the antibody (which was designed and commissioned by M.A.C. while at ProteoSys AG, [9,82]), corresponded to a CK2 consensus site. CK2 is a promiscuous kinase which gives rise to 10-20% of the human phosphoproteome [83]. On Western blots of breast cancer whole cell proteins, that antibody generated a low background smear of signals across many molecular weights, in addition to recognizing exogenously expressed PGRMC1 but not a S181A mutant (not shown). More rigorous characterization of the anti-pS181 antibody specificity would be required to exclude that the immunohistochemistry of tissue biopsies was not detecting increased levels of non-PGRMC1 CK2 phosphorylation sites.

Conclusion
Our work was inspired by the discovery that PGRMC1 phosphorylation differs between breast cancer subtypes differing in estrogen receptor expression, that also exhibit starkly contrasting levels of patient survival [9]. It is wholly feasible that characterization of the signal network surrounding PGRMC1 could lead to novel treatments for cancer, and other diseases involving PGRMC1 [2]. This work suggests that PGRMC1 is a key yet hitherto unrecognized protein involved in the early evolution of animals, and that its phosphorylation status is capable of influencing biological processes which can destabilize the ordered harmony of eumetazoan multicellular organismal collections of functionally specialized cells, with imperatively relevant connotations for human disease.

Establishment of cell lines and cell culture
Construction of cell lines and culture conditions is described in the companion publication [10]. Briefly, pcDNA3.1 plasmids expressing either a wild-type PGRMC1 open reading frame with a 3x heamaglutinin (3HA) tag (WT), an S57A/S181A double mutant (DM) or an S57A/ Y180F/S181A triple mutant (TM) were transfected into MIA PaCa-2 cells and stably selected by hygromycin selection. Three independent cell lines were established for each plasmid, and equivalent PGRMC1-HA expression levels for each line were established by Western blot. The exogenous PGRMC1-HA protein is present at approximately the same levels as endogenous MP cell PGRMC1. Individual cell lines have not been clonally selected, and so probably represent mixed lineages.

NpFR1 measurement of cytoplasmic redox status
Intracellular redox state measurement by fluorescent redox sensor NpFR1 was performed by flow cytometry as described for NpFR2 [10].

Genomic DNA isolation and sequencing
Genomic DNA for sequencing was isolated from one subline of each cell condition. Passaging was repeated every 4 days for thirty passages, with media changes every 48 h. To extract genomic DNA after thirty passages, approximately 5 × 10 6

MTT cell viability assays for effects of P4 on the dox and AG-205 responses
To assay for proliferative and protective effects of P4, 10 4 cells were seeded per 96 well plate well. For each cell type, n = 6 replicates were allowed to adhere for 3 h after which viable adherent cells were quantified by 3-(4,5-dimethylthiazolyl-2)-2,5-diphenyltetrazolium bromide (MTT) assay. Previous experiments had shown adherence to be complete after 1.5 h (data not shown). Identical replicates were incubated overnight in complete DMEM. After 23 h the media was exchanged for complete DMEM containing either 1 μM P4 or DMSO vehicle control. After a further 1 h incubation (t0) the media was exchanged for media containing either 1 μM P4 or DMSO vehicle control (to assay for proliferative effects, Fig. 2a), or with varying concentrations of doxorubicin (dox) (Fig. 2b-c) or AG-205, +/− P4 [36]. Cells were incubated for a further 24 h followed by viable cell quantification by MTT assay. The media was discarded, cells were washed with warm PBS and incubated with 100 μL of 0.5 mg/mL MTT (Sigma-Aldrich M2128) in phenol red free media (Sigma-Aldrich D1145). After 3 h incubation at 37°C and 5% CO 2 , the media was removed and 100 μL of DMSO (EMD Millipore, 317,275) was added to each well to solubilize the Formazan crystals. The cells were incubated for 1 h followed by mixing, and absorbance was read at 570 nm using a plate reader (Bio-Strategy P/L, Campbellfield, Vic., Australia). The percentage of viable cells was estimated by normalizing the absorbance of the treated or untreated cells to the average values from panel A.

Western blot
Cells were lysed with radioimmunoprecipitation assay buffer (RIPA buffer) (Sigma Aldrich, R0278) supplemented with protease and phosphatase inhibitor cocktail (Thermoscientific, 88,668). After scraping the cells, they were centrifuged at 8000 rpm for 20 min (Hermle Centrifuge Z233 M-2) at 4°C. 20 μg protein of thirty passage number of MP, WT, DM, and TM were loaded into the wells of a 12.5% SDS-PAGE gel and set to run at 150 V for 45 min on Power Pac. Transferring the protein from the gel to a PVDF membrane occurred under wet transfer with 1x Towbin buffer at 40 V for 3 h in an ice bath. The membranes were incubated with 1:2000 diluted Beta-Actin (Sigma Aldrich, AS441) and 1:2000 HA-Tag (Sigma Aldrich, H3663) primary antibodies in blocking buffer overnight at 4°C with shaking. Next day, the membranes were washed twice with 0.05% PBS-T and incubated with antimouse secondary antibody; horse radish peroxidase (HRP), for 1 h at room temperature. Next, the membranes were washed twice with 0.05% PBS-T and once with PBS. Colorimetric detection of the bands occurred using tetramethylbenzidine as described [10].

Hyperspectral fluorescence microscopy
UV and visible continuous wave epifluorescence microscopy with multiple excitation wavelength ranges from 335 nm to 532 nm and measuring emission in three emission wavelength ranges 450+/− 30 nm, 587+/− 17.5 nm and 700 nm long pass were used. Excitation wavelengths are supplied by low cost multi-LED light source optical fibre coupled to the Olympus IX71 microscope. An Andor/Oxford Instruments iXon 885 Electron Multiplying Charged Coupled Device (EMCCD) was used to capture images typically using an EM gain sufficient to lift the low light auto-fluorescence signal above the 17 electrons per pixel per second readout noise but with minimal contribution from clock induced charge. Gain linearity is ensured by using the Real Gain™ technology (Andor/Oxford Instruments). The camera is operated below − 70°C such that thermal noise is negligible. Quantum efficiency of the sensor varies from 50 to 65% across the range of interest from 450 to 670 nm and the intensity is digitized into~16 K values. Background signal is subtracted from all images which is kept minimal through the use of low fluorescence petri dishes (CELLview, Greiner Bio-One). A single wavelength image may take between 1 and 5 s depending on the sample and wavelength but 1-2 min is typically required for the entire stack of images for all spectral channels,, where the features for Fig. 1a were as follows.  Table S1. During growth under standard mid-log phase growth conditions, hyperspectral imaging was performed on~300 live cells per cell condition (MP, WT, DM, TM) and the mean cellular intensity in each channel was calculated. Pixel correlations between spectral channels were removed using principal component analysis (PCA), followed by a discriminatory projection to maximally separate the four cell groups (Fig. 1a). Technical details of this approach are described by Gosnell et al. [44]. In order to statistically analyze the cluster separations of Fig. 1a, an additional LDA projection of cell data was carried out, this time with two classes of cells chosen at a time. Their average spectra were projected onto a common direction identified by this additional LDA. The cell distributions were then tested by the Kolmogorov-Smirnov test.

Metabolite preparation
Scr-WT, scr-TM, shNNMT-WT, and shNNMT-TM cells as generated in "5.5. NNMT shRNA attenuation" above were seeded in Greiner Bio-One Tissue Culture Petri Dish 100 mm × 20 mm (Interpath, 664,160). The cells were rinsed with pre-warmed deionised water with brief shaking. The plates were placed on liquid nitrogen for 15 s. Cells were extracted with 1 mL of ice-cold extraction solvent (2:2:1methanol-ethanol-water) and harvested with scraper. The contents were transferred to 1.5 mL Eppendorf tube. The tubes were centrifuged at 4°C for 5 min at 16,100 x g and supernatant were transferred to a new 1.5 mL Eppendorf tube. The supernatant was filtered and evaporated dry under a gentle stream of nitrogen gas. The residual precipitates were resuspended in 500 μL of 60% acetonitrile/water (ACN/H 2 O) containing m-tyrosine (10 μg/mL) (Sigma, Castle Hill, NSW) as an internal standard.

Metabolite quantification
Metabolomic analysis was performed using an Agilent 1290 Infinity HPLC system equipped with a quaternary pump, degasser, temperature-controlled column and sample compartment coupled to an Agilent 6470 triple quadrupole (QQQ) mass spectrometer with a Jet Stream electrospray ionization (ESI) source (Agilent Technologies, Australia). The column temperature was maintained at 35°C and the autosampler temperature was set at 4°C. The electrospray (ESI) source settings were as follows: nebulizer gas, 45 psig; drying gas flow rate, 8.0 L/min; drying gas temperature, 300°C; sheath gas temperature, 350°C, sheath gas flow, 10 L/min; capillary voltage, 4000 V; and nozzle voltage, 0 V. The data were acquired positive ionization mode.
Methylnicotinamide samples and samples from Fig. 6a were analyzed with a Kinetex HILIC column (50 mm × 2.1 mm, 2.6 μm particle size, pore size of 100 Å, Phenomenex, CA, USA). The mobile phases were 100% HPLC-grade acetonitrile (A) and 10 mM ammonium acetate in HPLC-grade water adjusted to pH 3 with formic acid (B) at a flow rate of 200 μL/min. The solvent gradient began with 95% ACN for 0.5 min, decreasing linearly to 35% ACN over 12 min. The gradient was maintained at 35% ACN for 0.5 min and increased to 95% ACN over 0.1 min. The mobile phase composition was then held at 95% ACN for 10 min to re-equilibrate the column before the next injection. The total time for the gradient program was 24 min. Cellular extracts except Fig. 6a were analyzed with a HILIC-Z column (50 mm × 2.1 mm, 2.6 μm particle size, 100 Å pore size, Agilent Technologies, CA, USA) equipped with a guard column of the same stationary phase. The mobile phases were HPLC-grade acetonitrile:water (95:5) (A) and water (B), both containing 20 mM ammonium formate and 0.1% formic acid (pH~3) at a flow rate of 400 μL/min. The solvent gradient started at 100% A and was held for 0.5 min, followed by a linear gradient from 100% A to 60% A over 9 min. The gradient was maintained at 60% A for 0.5 min and returned to 100% A over 0.1 min. The solvent was then held at 100% A for 5.4 min to re-equilibrate the column before the next injection. Injection volume was 2 μL. The total time for the gradient program was 15 min.
Analysis of the set of samples was preceded and followed by injection of five concentrations of 1-MNA (10 μg/mL to 1 ng/mL in log steps) and six concentrations of S-adenosyl methionine (SAM) (100 μg/mL to 1 ng/mL) (Sigma-Aldrich, Castle Hill, NSW). In addition, after every 15 samples a blank injection and a mixture of 1-MNA, m-tyrosine and SAM was run to verify that there was no carryover between samples and that retention time and detector sensitivity had not changed during the course of running the samples.
The selective detection and quantitation of 1-MNA was achieved by monitoring the transition of m/z 137 → 108 using a fragmentor voltage of 120 V and collision energy of 20 V, with the further requirement that retention time matched that of the authentic compound (ca. 2.7 min). mtyrosine was quantitated by monitoring the transition m/z 182 -> 136 (fragmentor voltage 135 V, collision energy 20 V, RT = 4.9 min) and SAM was quantitated by monitoring the transition m/z 399 -> 136 (fragmentor voltage 100 V, collision energy 20 V, RT = 9.1 min). The fragmentor voltage and collision energy values were determined by optimizing system response to authentic 1-MNA, m-tyrosine and SAM. Limits of quantitation were determined by regressing the log(detector response) to log(concentration) response over the range of serial dilutions for each standard, and by observing where the dose/response curve departed from linearity.

Effect of 1-MNA on morphology
MiaPaCa-2 cells expressing DM PGRMC1-HA were seeded at 25% confluency and cultured overnight. The media was removed and the cells were incubated in the presence of the indicated amount of 1-MNA (Sigma-Aldrich, Castle Hill, NSW) for 24 h. The 1-MNA stock was made in sterile H 2 O. Five random images were taken per treated culture and assigned random numbers 1-20. These were given to a blinded scorer who decided for each image whether given cells were round or elongated.

Methylome assay
Genomic DNA from each cell subline were processed using the Illumina Infinium HD Methylation Assay (EPIC array) which interrogates > 850,000 CpG sites. The four cell lines were processed in triplicate equating to 12 samples, being three replicates of MP cells, and independent cell lines 1-3 for each of the PGRMC1 states [10]. The array was prepared at the Australian Genome Research Facility (AGRF) following the manufacturer's instructions. Quality checking of the samples was performed by Nanodrop Spectrophotometer and resolution on a 0.8% agarose gel at 130 V for 60 min. 500 ng total DNA was bisulfite converted with Zymo EZ DNA Methylation kit (Zymo Research, Irvine, CA, USA), following the manufacturer's standard protocol. Amplified DNA samples were fragmented and resuspended DNA was loaded onto a BeadChip. The BeadChip was incubated overnight while DNA fragments anneal. After hybridization the array was imaged on the Illumina iScan system.
Methylation data were processed using in R version 3.5.1 (www.r-project.org). Data were assessed using the hg19 build of the human referene genome annotated in the IlluminaHumanMethylationEPICanno.ilm10b4.hg19 bioconductor annotation. Quality control of the array was assessed with the R package Lumi [85]. Probes on the sex chromosomes, cross-reactive probes and probes with known SNPs at the CpG site were excluded. Following filtering 812,817 probes remained. The cleaned data were normalized with the SWAN normalization, and sample relationships examined. Genome-wide differential methylation analysis was undertaken using Limma [86]. Multiple testing correction was applied using the Benjamini and Hochberg method. Differential methylation in each group was observed on volcano plots and heatmaps.
Due to the excessive quantity of differentially methylated probes with small effects we refined the data used for functional analysis to the top 2000 most variable probes prior to pathways analysis. These probes explained > 90% of the total variance among the first two principle components while preserving a reasonable number of intragenic probes. Pathways enrichment was undertaken in R using the enrichPathway function in the ReactomePA package [87], which implements a onetailed hypergeometric test for overrepresentation. Multiple testing adjustment was applied with the Bonferroni & Hochberg FDR. KEGG pathways analysis was undertaken using the kegga function in the R package in LIMMA. Enrichment was calculated for genes from the top 2000 most variable probes with a p-value < 0.05. The gene sets were split to those with hypomethylated probes (logFC< 0) and those with hypermethylated probes (logFC> 0) for enrichment. Enrichment of Gene Ontology (GO) terms was undertaken using the GOstats package [88]. Ontology assignments are calculated using a hypergeometric test for overrepresentation on the significantly differentially methylated probes (p < 0.05) among the 2000 most variable. Multiple testing adjustments were applied using the FDR. Intersects of common pathways from KEGG and Reactome analyses were assessed using the UpSetR package [89].
after treatment by an NNMT-specific shRNA (NN) or a random scramble shRNA control (scr) in TM cells stably expressing lentiplasmid-driven shRNAs. RT-PCR Methods follow A. p < 0.003 (2-tailed T test). Methods follow was applied separately to significant hypermethylated and hypomethylated probes. The main bar chart presents the number of pathways common to each comparison. Each bar represents an intersect notated by linked dots on the x-axis. The number of significant pathways per comparison is presented as a horizontal bar chart. The KEGG results are in File S1 and the list of pathways unique to each comparison are provided. Pathways analysis of probes genes corresponding to the top 2000 most variable differentially methylated probes in CpG Island or Shore chromosomal regions. Probes were separated into hypermethylated and hypomethylated data sets from D, including only probes from Island or Shores. Enriched Gene Ontology (GO) pathways enrichments are shown. The identities of pathways uniquely differential to one of the three cell type comparisons (WT/MP, DM/WT, or TM/DM) for both hypoand hyper-methylated data sets and KEGG pathways are available in Supplemental Information File 1.
Additional file 8 Table S1. Spectral channel filters employed for hyperspectral autofluorescence imaging. Related to Fig. 1. For further details of this approach see Gosnell et al. [44].
Additional file 9 Table S2. GO pathways enrichment results. Related to Fig. 7. The top ten GO enrichments for each cell comparison from Fig. 7. Full results are available in File S1.
Additional file 10 Table S3. Unique KEGG pathways detected for TM/ DM hyper and hypo-methylated gene sets. Related to Fig. 7. The identities of the 11 hypermethylated and 7 hypomethylated pathways unique to the TM/DM comparison from Fig. S7B are given. Full results are available in File S1.
Additional file 11 File S1. A zip archive containing excel files with the most significant pathways enrichment results from methylomics analysis including those of Fig. S7B, Table S2 and Table S3. The archive unpacks as five separate folders, each containing the following pathways analysis results. A) All GO enrichments for each cell comparison significant below the adjP =0.001 level. B) All GO enrichments that were uique to a particular cell comparison. C) KEGG enriched pathways for each cell comparison, with analyses performed separately on hypermethyalted (files labelled up) and hyperemethylated (files labelled down) probes. These correspond to the pathways of of Fig. S7B. D) All KEGG enriched pathways from C which were unique to a particular cell comparison. E) All All Reactome enriched pathways for each cell comparison significant below the adjP =0.05 level.