Transcriptomics uncovers substantial variability associated with alterations in manufacturing processes of macrophage cell therapy products

Gene expression plasticity is central for macrophages’ timely responses to cues from the microenvironment permitting phenotypic adaptation from pro-inflammatory (M1) to wound healing and tissue-regenerative (M2, with several subclasses). Regulatory macrophages are a distinct macrophage type, possessing immunoregulatory, anti-inflammatory, and angiogenic properties. Due to these features, regulatory macrophages are considered as a potential cell therapy product to treat clinical conditions, e.g., non-healing diabetic foot ulcers. In this study we characterized two differently manufactured clinically relevant regulatory macrophages, programmable cells of monocytic origin and comparator macrophages (M1, M2a and M0) using flow-cytometry, RT-qPCR, phagocytosis and secretome measurements, and RNA-Seq. We demonstrate that conventional phenotyping had a limited potential to discriminate different types of macrophages which was ameliorated when global transcriptome characterization by RNA-Seq was employed. Using this approach we confirmed that macrophage manufacturing processes can result in a highly reproducible cell phenotype. At the same time, minor changes introduced in manufacturing resulted in phenotypically and functionally distinct regulatory macrophage types. Additionally, we have identified a novel constellation of process specific biomarkers, which will support further clinical product development.


Macrophage differentiation methods. As schematically presented in
, purified monocytes were differentiated into M0, M1, M2a, Mreg type-of-cells or PCMO-like cells according to the published protocols with minor modifications.
Mreg cells were differentiated from purified monocytes for 7 days similarly to Mreg_UKR cells, but in Cell + flasks (Sarstedt) 1,30 and with the following modifications: (1) medium was changed twice during the process (day 1 and 4); and (2) cells were harvested by scraping.
PCMO-like cells were cultured and harvested similarly to Mreg cells, but in the presence of 0.4 ng/ml of recombinant human IL-3 (Miltenyi Biotech) and without addition of IFN-ɣ on day 6. This production method deviates from the one originally reported for PCMO 19,32 to better fit the scheme and timelines of the manufacturing of other cells. Because of these changes in manufacturing, we call these cells PCMO-like throughout the text.
M0, M1 and M2a cells were produced according to protocols given in Hutchinson et al. 30 in the gas-permeable MACS GMP bags. Cells were differentiated for 6 days in RPMI 1,640/GlutaMAX/antibiotics (same as for Mreg_ UKR), containing 5 ng/ml M-CSF and 20% foetal bovine serum (FBS; Gibco). On day 6 medium was exchanged, serum concentration reduced to 5% FBS and M-CSF increased to 25 ng/ml. M1 cells were polarized with 100 ng/ ml of lipopolysaccharides (LPS, Sigma-Aldrich) and 25 ng/ml IFN-γ, while M2a cells were polarized with 20 ng/ ml IL-4 (R&D Systems). M0 cells did not receive any additional cytokines. Cells were harvested on day 7. 200 system (Bio-Rad). Data were analyzed with Bio-Plex Manager software version 4.1.1, and visualized utilizing R. In order to compare multiple assays, median fluorescence intensities (MFI) were first normalized (nFI) by subtracting the background measurement from each analyte. Principle coordinates were calculated in R in order to assess the degree of inter-assay technical variability in comparison with experimental, inter-cell-type variability. Due to the inter-assay variability, a reduced resolution, binned heatmap method was employed. The R package dr4pl v1.1.11 34 was used to create a combined standard curve from which the binning thresholds were established. Bins were drawn from the upper and lower asymptotes, as well as the calculated midpoint (EC50).
Gene expression characterization by RT-qPCR. Total RNA was extracted from cells using RNeasy Protect Cell Mini Kit (Qiagen) with QIAshredder homogenizers (Qiagen). RNA was treated with DNase I (PrimerDesign or Invitrogen) according to manufacturer's instruction, samples' concentration was adjusted to 50 ng/µl and verified on Qubit Fluorometer (Invitrogen) using Qubit RNA HS Assay Kit (Invitrogen). Expression of DHRS9 and Ido1 mRNA was measured by RT-qPCR relative to GAPDH mRNA endogenous control using TaqPath 1-Step Multiplex Master Mix Kit (Applied Biosystems) and QuantStudio 5 Real-time PCR system (Applied Biosystems) according to manufacturer's instruction. TaqMan assays were from Applied Biosystems: for Human Ido-1 assay number Hs00984148_m1 (FAM-MGB), for Human DHRS9 Hs00608375_m1 (FAM-MGB) and for Human GAPDH Hs03929097_g1 (VIC-MGB). Each amplification reaction contained 50 ng of RNA and was performed in triplicates. Data was analyzed with QuantStudio Design and Analysis desktop Software (Applied Biosystems). Changes in expression (Rq) were calculated relative to expression of corresponding mRNA in the starting material, i.e. CD14 + monocytes. phagocytosis assay. The flow cytometry-based phagocytosis assay was performed with harvested macrophages, CD14 + monocytes and CD3 + T cells using the pHrodo Green E. coli BioParticles Phagocytosis Kit for Flow Cytometry (Invitrogen) according to the manufacturer's protocol. One million macrophages/monocytes and half a million of T cells were incubated with pHrodo Green E. coli BioParticles Conjugate for 1 h at + 37 °C or on ice in order to inhibit particle uptake. Cells were harvested and stained with 7-AAD, CD45-PE/Cy7 (Biolegend), and CD3-BV421 (Biolegend) for T-cell samples. For analysis, a minimum of 4 × 10 4 live cells (doubletdiscriminated CD45-positive, 7-AAD-negative) were acquired and analyzed in the FlowJo software (TreeStar). Phagocytosis was determined as an increase in pHrodo Green fluorescence relative to the one observed with control sample, incubated on ice, and calculated as the percentage of pHrodo Green-positive cells.

RNA-Seq sample preparation. Enriched monocytes and macrophages were harvested directly into
RNAprotect Cell reagent (Qiagen) and frozen at − 80 °C. RNA extraction and sequencing were performed by Eurofins Genomics to a minimum of 30 M reads using Illumina HiSeq 4000 with 2 × 150 bp paired-end sequencing procedure.
RNA-Seq data analysis. Quality control. Quality of the RNA-Seq reads was inspected using MultiQC v1.7 35 .
Read alignment and read counts. RNA-Seq reads were aligned the using STAR aligner v2.5.2 36 , against the Ensembl GRCh38 patch 12 release 95 of the human genome. The aligned reads were counted using the Bioconductor SummarizeOverlaps in 'union' mode. Reads were tallied from the reverse strand, using the hidden argument 'preprocess.reads = invertStrand' . Under-sequenced samples were combined post-quantitation, using the DESeq2 collapseReplicates function of DESeq2 R package 39 v1.24.0. Read counts were normalized using the variance stabilizing transformation (VST).
Hierarchical clustering analysis. Ensembl gene VST values from the sequenced samples were clustered by Ward2 clustering using the Pearson correlation coefficient as a distance metric. Trees were resampled to 1,000 bootstrap iterations using pvclust v2.0-0 40 .
Visualisations of genes. Volcano plots of differentially expressed genes were created using ggplot2 43 . Gene expression heatmaps display VST-transformed values and are plotted using pheatmap v1.0.12.
Enrichment analysis. Gene set enrichment analysis was performed using the gene sets as collected under MSigDB v7.0 44,45 . Enrichment was calculated using fgsea v1.10.1 46 . The Reactome 47 and the immunologic signatures 48 datasets were selected for presentation. Enrichments with |NES|> 1 and adjusted p-value < 0.01 were considered significant. The Database for Annotation, Visualization and Integrated Discovery (DAVID 6.8) 49,50 was also used to identify overrepresented biological terms in KEGG pathways and Gene Ontology 51,52 set enrichments. Protein-protein interaction networks in specific gene sets were analyzed using STRING 53 .
Further gene set enrichment, networks, predicted activations, and functional analyses were performed using Qiagen Ingenuity Pathway Analysis (IPA) 54  www.nature.com/scientificreports/ ggplot2. "Canonical pathway" activation predictions were considered significant if they exhibited a p-value less than 10 -4 . "Biofunction activation" predictions were considered significant if they exhibited a B-H corrected p-value < 0.05 and a bias-corrected z-score > 2. Predictions that referred to a disease, cancer, abnormality, or injury state were ignored. Results which were flagged for bias had their bias-corrected z-score substituted for the uncorrected value, per vendor documentation. Statistical analysis. Statistical analyses for assays other than RNA-Seq were performed using GraphPad Prism 5.04. Results are expressed as mean ± SD. Comparisons between groups were performed using one-way ANOVA with Tukey's or Dunnett's post hoc test. P < 0.05 was defined as statistically significant. The relationship between fold-changes in gene expression in RNA-Seq and RT-qPCR analysis was investigated using linear regression.
ethics. Peripheral blood monocytes were isolated from healthy donors' leukapheresis products. Informed consent was obtained from all donors in accordance with the Declaration of Helsinki. The studies were approved by the Research Ethics Committee of the Northern Savo Hospital District (approval number 186/2017).

Results
Generation and initial characterization of macrophages. The three potential cell therapy products, Mreg, Mreg_UKR, and PCMO as well as comparator macrophages, M0, M1, and M2a were differentiated from donor-derived monocytes according to published protocols as outlined in Fig. 1a. All products were assessed by flow cytometry for 23 surface markers. Only a few of the screened extracellular proteins were suitable to specifically identify each cell type (Supplementary Figure S1). Therefore, we utilized a principle component analysis (PCA) of Wisconsin-standardized nFI consisting of all the 23 markers, and identified each cell type with better confidence (Fig. 1b).
The most prominent differences between each cell type were observed within the levels of activation-associated, pro-inflammatory markers. In accordance with the subsequent RNA-level findings (see below), CD38 and CD40 were subtly upregulated in Mreg cells, while M1 macrophages displayed markedly elevated levels of CD38, CD40 and CD80 (p < 0.001 for CD38 and CD80, and p < 0.05 for CD40), when compared to other macrophage types (Supplementary Figure S1). However, this limited set of functionally related markers aside, only a few others, e.g. CD71 for Mreg_UKR, offered a reliable way to identify a specific cell type if employed alone. This is not surprising due to the shared starting material, as well as the known plasticity and overlapping functionality of different macrophage phenotypes 27,36,37 . When used in conjunction with one another the full panel of markers provided sufficient information to unsupervised machine learning techniques that linear discriminant analysis (LDA) (Fig. 1c) and Wisconsin-standardized PCA were able to discern between most cell types (Fig. 1b). Only the M0 cell type defied attempts at classification. The resultant ordination shows each cell type cluster arranged radially around a highly variable M0 core, with each component loading vector, indicating the relative importance of each marker in informing the location of each sample in a cell type cluster.

Unbiased classification of macrophages by RNA-Seq.
For an in-depth, unbiased characterization of the produced macrophages, the transcriptome of the cells was assessed by RNA-Seq. As expected, all macrophage types were clearly distinct from monocytes by principal component analysis (PCA) and hierarchical clustering of the samples (Fig. 2). The first principal component (PC1 in Fig. 2a) explained 65% of the variance mainly driven by the difference between monocytes and all macrophages that were produced. The PC2 separated the differentiated macrophages and explained 18% of the variance, with M1 (pro-inflammatory) and M2a (wound healing) macrophages showing the furthest separation. Consistent with PCA analysis, the highest number of differentially expressed genes (DEGs), > 5,000, was observed between monocytes and each of the produced macrophages, followed by the M1 vs M2a macrophages with over 3,000 DEGs (Table 1). Together, Mreg_UKR and Mreg formed a common cluster of regulatory macrophages between M1 and M2a. The PCMO-like cells were highly similar to M2a and M0 macrophages and dissimilar to the regulatory macrophages as reflected in PCA and in the low number of DEGs observed between M2a and PCMO-like cells vs non-polarized M0.
Exclusion of CD14 + cells in the gene expression PCA plots facilitated a better separation of all the macrophage types (Fig. 2b). The variation between cell types was always greater than the variation between the donors in the same cell type cluster, indicating the robustness and reproducibility of the manufacturing methods (Figs. 1a and 2b). Although the two regulatory macrophage products, Mreg and Mreg_UKR, were most related to each other, they still had clearly distinct transcriptomes with 1,002 DEGs between them. Nevertheless, the differentiating medium and added cytokine (IFN-γ) were identical. The difference between Mreg and Mreg_UKR appeared to be solely driven by the culture vessel and medium changes in the Mreg process. The hard surface of the flasks and medium replenishment made Mreg cells somewhat more M1-like, when compared to Mreg_UKR. Accordingly, in flow cytometry analysis, slightly elevated levels of CD38 and CD40 were noted for Mreg cells (Supplementary Figure S1). Additionally, just as M1 macrophages, Mreg cells secreted chemoattractants MCP-1, MIP-1a and MIP-1b, which play a role in recruiting other cells to sites of inflammation, but not classical inflammatory cytokines like TNF-α, IL-6, or IL-12 nor pro-inflammatory chemokine RANTES, which were secreted only by M1 (Supplementary Figure S2). From regulatory macrophages, Mreg_UKR were more distinct from M1 and M2a macrophages than Mreg as reflected in the total amount of observed DEGs and protein-coding genes (Table 1).
Hierarchical clustering confirmed that all samples, while distinct, were most closely associated with samples of the same cell type (Fig. 2c). M2a, PCMO and M0 formed a separate branch from cells that were produced with Scientific RepoRtS | (2020) 10:14049 | https://doi.org/10.1038/s41598-020-70967-2 www.nature.com/scientificreports/ the addition of IFN-γ: M1, Mreg and Mreg_UKR. Within the IFN-γ-stimulated branch, regulatory macrophages clustered separately from M1 cells. Taken together, the consistency of the manufacturing process and the reproducibility of manufactured cells' phenotypes were ascertained by the hierarchical clustering and PCA ordination.

Characterization of comparator phenotype macrophages and PCMO-like cells. The transcrip-
tomes of M1, M2a and PCMO-like cells were interrogated for expression of their literature-reported phenotype markers. Furthermore, the M1-and M2a-polarized macrophages were compared to each other, and to the nonpolarized M0 to identify additional genes induced in response to polarizing agents, IFN-γ/LPS (M1) and IL-4 (M2a). The plethora of previously reported polarization markers was narrowed down to a set of the most widely accepted human marker genes 38 and complemented with the genes that were reported in previous M1/M2a comparative transcriptome analyses 9,23,25,27 . The full list of examined genes and their expression in produced M1 (162 genes) and M2a (123 genes) macrophages is given in Supplementary Table S2. Overall, the transcriptional profiles of M1 and M2a macrophages conformed to those described in the literature ( Fig. 3a and Supplementary  Table S2).
To screen for novel M1/M2a polarization markers, we analysed the genes with 16-fold upregulation (LFC4) in M1 and M2a relative to each other and to M0 macrophages. This lead to identification of 224 protein-coding transcripts highly expressed in M1 relative to both M2a and M0 and of 15 highly abundant transcripts in M2a relative to both M1 and M0 (Fig. 3b-e). Of interest were 25 transcripts highly upregulated in M2a in comparison to M0 macrophages (M2a_vs_M0) (Fig. 3b, Table 2). Of these transcripts, 14 were not on our compiled marker  Table S2). Analysis of functional interactions between the unreported M1 145 marker genes ( Table 2) using STRING revealed a tight gene network centred around TLR7/JAK3/IFIT1 genes which are known activators [39][40][41] and regulators 42 of pro-inflammatory response in macrophages with enrichment for the terms: "immune response" (45 genes to GO:0006955), "response to LPS" (14 genes to GO:0032496) and "immune system process" (51 genes to GO:0002376) (Fig. 3e). Unlike with M1 and M2a macrophages, we found little resemblance of the manufactured PCMO-like cells to the PCMO described in literature. The PCMO-like cells were produced as described previously for PCMO 19,32 apart from using a highly purified CD14 + population and culturing for a longer time period (7 instead of 4-6 days) without β-mercaptoethanol. Possibly, because of these modifications, no specific reactivation of pluripotency genes (POU5F1 (OCT4), NANOG and MYC) 18 , nor elevated expression of THY1 (CD90), CSF1R (CD115) and IL3RA (CD123) 43 were observed when compared to either monocytes or to other produced macrophages (Supplementary Table S3

Regulatory macrophages show intermediate phenotype between M1 and M2a. M1 or M2a
polarization resulted in a massive shift in the global gene transcription pattern when compared to non-activated M0 macrophages. Perturbation involved significantly more genes for M1 (3,281) than for the M2a (621). It appears that the default macrophage differentiation, in M-CSF containing medium 9,36 , is driven more towards an M2a phenotype even without the presence of polarizing cytokines. Therefore, comparing M2a vs M1 macrophages in relation to M0 partially masks drastic differences between M2a and M1 transcriptomes (Fig. 3a).
Because of that, direct a M1_vs_M2a comparison was used to analyse how regulatory macrophages or PCMOlike cells fit into the M1-M2a polarization spectrum on the individual gene expression level. Initially, the 46 most upregulated (M1-specific, basemean > 1,000) and the 47 most downregulated (M2a-specific, basemean > 200) genes were selected and their expression was analysed in regulatory macrophages and PCMO-like cells in relation to M1 and M2a by heatmaps (Fig. 4a). IFN-γ stimulated genes were also induced in Mreg and Mreg_UKR, but to a lesser extent than in M1 macrophages. Nevertheless, some IFN-γ-induced genes, e.g. CXCL10, CXCL11 and IL32, expected to be upregulated in Mreg_UKR remained at a similar expression level to M2a and PCMOlike cells. On the other hand Mreg, Mreg_UKR and PCMO-like cells shared increased expression of some www.nature.com/scientificreports/ M2a-specific genes (Fig. 4a). The expression of mannose receptor MRC1 (CD206), which is one of the key markers of the M2a phenotype, was higher in Mreg_UKR than in M2a and further supported by flow cytometry (Supplementary Figure S1). While some IL-4-induced genes, MRC1, PPARG and FABP4, were also upregulated in Mreg_UKR, Mreg, and PCMO-like cells, the expression of the others, e.g. ALOX15 and FCER2, was induced only in M2a. The increased expression of the first subset was also noted in unstimulated M0, meaning that the default macrophage maturation likely drives expression of those genes, which are then turned off in the presence of the strong M1 polarization signal such as IFN-γ plus LPS, but not IFN-γ alone. Next, the DEG data sets were analysed by Ingenuity Pathway Analyser (IPA) to assess changes in the most significant molecular networks, biological functions and canonical pathways between all cell types pairwise. IPA upstream regulator analysis of M2a_vs_M1 and M1_vs_M2a comparisons uncovered key transcriptional  Table S5). Visualization of the IPA-generated gene networks driven by the identified M2a upstream regulators revealed high expression of the target genes in M2a when compared to M1 macrophages (Supplementary Figure S4a). Similarly, genes in the networks for the identified M1 regulators were highly expressed in M1 when compared to M2a (Supplementary Figure S4b). Both M1 and M2a regulator gene networks were only partially upregulated in Mreg and Mreg_UKR, with a subset of genes in M1 regulator network downregulated in Mreg_UKR (Supplementary Figure S4a,b). The PCMO-like cells showed high expression for a great proportion of the IL-4, IL-13 and TGFB1 regulated genes, but displayed either downregulation or no differential gene expression relative to M2a for M1 regulatory networks. Therefore, major network analysis also supported the intermediate phenotype for Mreg and Mreg_UKR between M1 and M2a, while PCMO-like cells were most similar to M2a.
Anti-inflammatory and antimicrobial biofunctions are predicted for regulatory macrophages. Additionally, using IPA, we examined the main perturbed pathways and biofunction-activation scores in comparison with M0, used here as a common baseline phenotype. The IPA predicts activated and inhibited biofunctions based on z-score value. We selected activated and inhibited biofunctions with │z-score│ ≥ 2 for M1 and M2a macrophages and investigated predictions for the same biofunctions in Mreg, Mreg_UKR and PCMO-like cells (Fig. 4b). The majority of the predicted activated biofunctions for M1 were related to cell activation, cell migration, inflammatory and antiviral/antimicrobial response. For M2a the main activated biofunctions were related to cell migration. According to IPA analysis Mreg, Mreg_UKR and PCMO-like cells were predicted to behave differently on a number of biofunctions when compared to both M1 and M2a macrophages (Fig. 4b). Specifically, cell movement, migration and chemotaxis z-scores were negative for Mreg, Mreg_UKR and PCMO-like cells. Cell activation z-scores were strongly negative for Mreg_UKR and partially for PCMO-like cells as opposed to M1. Inflammatory response z-scores were negative for Mreg_UKR, Mreg and PCMO-like cells as opposed to M1 and exceeded that of M2a. The z-scores for antiviral and antimicrobial response for both Mreg and Mreg_UKR were in concordance with M1. Thus, based on global transcriptional profile Mreg_UKR, Mreg and PCMO-like cells were predicted to be functionally different from M2a and M1 macrophages. This is specifically intriguing for PCMO-like cells, which based on individual gene expression were similar to M2a and M0, but according to IPA were more anti-inflammatory, with less cell activation and migration properties than M2a and M0. At the same time, Mreg_UKR and Mreg were more anti-inflammatory than the comparator macrophages possessing antimicrobial properties and low cell activation and migration properties. Functional characterization of DEGs in IPA identified numerous, significantly enriched canonical pathways (FDR-adjusted P value ≤ 0.05): 222 pathways for M1_vs_M0, 41 for M2a_vs_M0, 135 for Mreg_vs_M0, 116 for Mreg_UKR_vs_M0, and 43 for PCMO_vs_M0. We selected the most significantly enriched pathways for M1 and M2a with an FDR-adjusted p < 0.00001 and compared the enrichment significance of the same pathways in Mreg, Mreg_UKR, and PCMO-like cells (Fig. 4c). The Antigen Presentation, Complement system, T helper cells differentiation and activation, and OX40 signalling pathways were more significantly enriched in regulatory macrophages than in M1, suggesting that they could be more effective in these functions.  Table S6 and Supplementary Figure S5a). Hierarchically clustered heatmap analysis for the literature-derived gene set revealed concordant gene expression in Mreg and Mreg_UKR for the majority of these genes. The notable exceptions were CCL1 38 , F13A1, XCR1, TNFSF18, LYVE1 and CX3CR1 2 , which were not expressed in our regulatory macrophages. In agreement with the literature, and confirmed by our flow cytometry analysis, expression of mRNA for co-stimulatory molecules CD80 and CD40 was lower in Mreg and Mreg_UKR than in M1. Mreg and Mreg_UKR were generally most positive for expression of immunoglobulin Fc receptor, FCGR1A, scavenger receptor MSR1 and DHRS9, which is a recently proposed marker of regulatory macrophages 44 . Nevertheless, with exception for FCGR1A, none of these markers exhibited a unique expression pattern that could be used to define regulatory macrophages exclusively and pass a statistical significance test. Indeed, when specifically  Figure S5b,c). We endeavoured to pinpoint prospective markers for regulatory macrophages from our data. Differential expression gene sets were analysed using Venn diagrams to identify genes specifically upregulated (LFC1) in both Mreg and Mreg_UKR in comparison to M1, M0, M2a and PCMO-like backgrounds. In both regulatory macrophages, 16 genes were commonly upregulated ≥ twofold, with a further 27 genes in Mreg and 60 genes in Mreg_UKR (Fig. 5a). Protein products from the majority of these genes are known to localize to the cell membrane, with specific enrichment for molecules involved in antigen presentation (p < 0.001). Hierarchically www.nature.com/scientificreports/ clustered heatmap analysis of expression levels for those genes revealed that only a subset was exclusively upregulated in either Mreg or Mreg_UKR (Fig. 5b-d).
Similarly, we defined genes specifically downregulated in regulatory macrophages when compared to other products. Only 4 genes, FOXP1, ZNF827, MAOA and XYLT1 were identified in Mreg when compared to M1, M2a, M0 and PCMO-like cells. MAOA and XYLT1 belong to the TGFB1 regulatory network and are associated with the M2 phenotype 27,45 . The corresponding list of Mreg_UKR-specific downregulated genes contained 122 protein coding genes. Protein products for the majority of those genes are annotated as being localized to the cell membrane (73 genes to cellular component GO:0016020 membrane), with 15 genes involved in inflammatory response (GO:0006954) and 19 genes involved in cell migration (GO:0016477) (Supplementary Figure S5d). Inflammation associated genes IL1B, CCL3 (MIP-1) and CCL4 (MIP-1b) were further verified by multiplex ELISA and shown to be secreted at negligible levels in Mreg_UKR (Supplementary Figure S2a). While the regulation of cytokine secretion is a spatiotemporal process, the cytokine panel was nevertheless transformed to log 2 -fold changes and compared to RNA-Seq estimates in order to investigate a possible correlation. Concordance was achieved for a majority of the panel (Supplementary Figure S2b), with only a small subset of cytokines, including IL-1β, not secreted proportionally to the number of expressed transcripts (Supplementary Figure S2c).
Based on these data, we have identified a number of genes specific for regulatory macrophages and further distinguishing the two studied regulatory macrophages, namely Mreg and Mreg_UKR.

Discussion
Cell therapy applications are promising treatments for a variety of diseases. However, success in the fundamental research is rewarded with new challenges while moving towards clinical trials. In cell therapy, the process makes the product and, consequently, the development of GMP-compliant manufacturing and supportive analytics must occur before pre-clinical and clinical studies. An ample analysis toolbox should be established based on an unbiased characterization of the developed product to secure process reproducibility and the resultant product's potency.
In this study, we sought to deepen the knowledge of regulatory macrophages, their polarization states and the impact of manufacturing on their phenotype. We adapted earlier published manufacturing methods for regulatory macrophages and PCMOs, and assessed their characteristics relative to the comparator macrophages (M0, M1, and M2a). We focused on a set of well documented macrophage markers determined by flow cytometry 1,2,31,48 . However, the produced cell types were not conclusively distinguished using these markers. Furthermore, we wished to understand the consequence of any process change in order to ensure that observed shifts in phenotype do not have an undesired effect on the functionality and/or safety profile of a cell therapy product.
The necessity of an unbiased survey of gene expression for cell therapy products has been emphasized in a study comparing bone marrow stromal cells (BMSCs) produced at different centres according to internal protocols 49 . Hierarchical clustering and principal component analysis of RNA-Seq data showed that variability in manufactured BMSC was greater between different centres than within each centre. This inter-centre discrepancy likely resulted in variation in the BMSCs functions, such as the therapy's ability to form bone and support haematopoiesis in an in vivo transplant model.
We employed RNA-Seq for comprehensive, unbiased molecular characterization of macrophages included in this study and were able to reliably distinguish all macrophage types. Importantly, RNA-Seq data confirmed that the manufacturing processes were highly reproducible and largely unaffected by donor-to-donor variations. The inflammatory M1 were the most distinct macrophage type that could be also distinguished utilizing few phenotype markers. The non-polarized M0, alternatively activated M2a, and PCMO-like cells had minor differences as visualized in principal component analysis and hierarchical clustering. M-CSF is common in manufacturing of all these cell types and appears to have a greater effect on the cell characteristics than does the source of serum or addition of IL-3 or IL-4. Clinical and in vivo studies suggest that PCMO, much like M2 macrophages, could facilitate tissue remodelling and wound healing 19,50 .
The two regulatory macrophage products clustered apart from each other based on RNA-Seq data analysis. The original regulatory macrophage cell therapy product was developed as an adherent Mreg culture in flasks 30 . Later, the production was up-scaled to a GMP-compliant process (Mreg_UKR), using semi-adherent bag conditions 31 , and in parallel, possibly to minimize handling steps, medium replenishments were abandoned. Macrophages are highly attuned to their microenvironment, and we show here that the putatively minor changes in the regulatory macrophage manufacturing resulted in a clear phenotypic shift which was reflected in surface protein expression, transcriptomics, attenuation of IFN-γ signaling pathways, and activation status reflected in secretome. Considering the plasticity of macrophages it remains to be investigated whether these in vitro differentiated regulatory macrophages will retain their characteristics in a complex in vivo environment.
IFN-γ is a strong activating cytokine which is used to polarize regulatory and M1 macrophages. It activates the expression of number of genes including indoleamine 2,3-dioxygenase (IDO1), IL32 and CXCL10 (IP-10). However, these genes and corresponding proteins were less expressed in Mreg_UKR than in Mreg, despite receiving an identical IFN-γ dose. Regulatory macrophages exert their immunomodulatory function by transforming naïve CD4 + T cells into anti-inflammatory regulatory T cells 2 and by suppressing T cell proliferation, both of which are dependent on increased IDO activity 1 . The proposed mechanisms include direct elimination of effector T cells by tryptophan deprivation 1,47 and indirect immunoregulation via Stat3-mediated induction of regulatory T cells 51 . We observed a clear difference in IDO1 expression between manufacturing conditions. In the Mreg_UKR, which are produced without medium changes, IDO1 expression is diminished on both mRNA and protein level. Although Mreg_UKR manufacturing resulted in somewhat "dormant" cells, they retain their propensity for re-activation in a different microenvironment. Our preliminary data suggest that upon stimulation with LPS Mreg_UKR not only upregulate expression of IDO1, but also secrete wound healing cytokine IL-10 (unpublished data). This further supports an anti-inflammatory and immunomodulatory role for Mreg_UKR, which in the chronic wound setting could contribute to resolution of inflammation rather than exacerbate the immunopathology.
A recent study demonstrated that regulatory macrophages not only secrete angiogenic proteins, but this secretion increases under induced hypoxia 52 . We did not detect any significant gene expression or secretion of cytokines/growth factors associated with angiogenesis in regulatory macrophages as compared to other phenotypes under basal conditions (i.e. right after harvest). Potentially minor differences in manufacturing or experimental design could account for the observed disagreements.
Determination of objective and measurable biomarkers that would reflect potency are in high demand for potential cell therapies 53 . This study provides a foundation to develop a toolbox for in vitro and in vivo phenotype and potency analyses. Although functional studies were limited to secretion and phagocytosis assays, the pathway and biofunction analysis correctly predicted the known functions of M1 and M2a macrophages. Functional predictions for regulatory macrophages were distinct from those for M1 or M2a and prognosticated Mreg and Mreg_UKR to have anti-inflammatory and antimicrobial attributes. Such properties would be essential in a macrophage cell therapy product that seeks to resolve inflammation in e.g. chronic wound, where the reduced Scientific RepoRtS | (2020) 10:14049 | https://doi.org/10.1038/s41598-020-70967-2 www.nature.com/scientificreports/ capability of macrophages to phagocytose debris and apoptotic neutrophils results in prolonged inflammation 54 . The membrane-localized proteins (identified here by RNA-Seq), especially the ones involved in antigen presentation and phagocytosis, HLA-DPA1, HLA-DPB1, HLA-DOA and CD64 (FCGR1A), are the prime choices for development of novel flow cytometry markers for fundamental research and release panels reflecting both phenotype and therapeutic potency of macrophage products. As part of a release panel, specific cytokines, such as IP-10 for Mreg, could be measured by ELISA directly from the culture medium without sacrificing the cell numbers in the final macrophage product. Another option is to develop a quantitative, multiplex RT-qPCR assay or targeted RNA sequencing panel for rapid analysis of several genes simultaneously. Such a panel could be based on the cell type specific transcripts identified in this study, e.g. the ones encoding membrane-localized proteins and specific transcription factors, such as SMARCD3, CEBPE and VENTX for Mreg_UKR. The role of these transcription factors in promoting the regulatory macrophage phenotype has yet to be investigated. Interestingly, VENTX is known to upregulate expression of a number of genes, including CD71, CD64, CD206 55 , which showed increased expression on flow cytometry and/or in RNA-seq. Selection of target transcripts can be further assisted by machine learning based on RNA-Seq data, which, as shown here, is a powerful and operatorindependent tool in macrophage classification.
To our knowledge, we provide the first comprehensive analysis of clinically relevant regulatory macrophages and report several unexpected differences. Our transcriptome analysis showed that this method not only distinguishes the different macrophages from each other, but also reveals new potential markers and predicted biofunctions that could be used in the determination of functionality, potency, and ultimately in a product release panel. Process reliability and reproducibility are the cornerstones of successful cell therapy, reducing the risk of changes in cell function which may compromise the safety or efficacy of the final product. We have demonstrated that the manufacturing processes were highly consistent and reproducible. Simultaneously, we show that even seemingly minor process changes affected the transcriptome of the regulatory macrophage product. Further studies could address whether these transcriptomic fluctuations translate into a functionally different cell therapy in vivo. Our results emphasize the importance of early process development for translation into clinical settings as well as the importance of concurrent unbiased product characterization to ensure reproducibility of the desired phenotype, functional state, and designation of final release markers for the end product.