Transcriptional Analysis Reveals Alteration of Lung Macrophages Following Mesenchymal Stem Cell Therapy in Neonatal Mice With Hyperoxia-Induced Injury


 Background Lung immaturity is one of the most serious consequences of growth restriction and premature birth. Preterm babies often require mechanical ventilation to survive, but exposure to high levels of oxygen can permanently damage the lungs and interrupts normal development. As lung macrophages play an important role in hyperoxic lung injury and repair, our objective was to use next generation sequencing (NGS) to identify changes in the macrophage transcriptome following neonatal hyperoxia, with and without treatment with human mesenchymal stem cells (hMSCs). We provide the first report of RNA-sequencing of lung macrophages following neonatal hyperoxia and hMSCs therapy. Methods Neonatal mice exposed to normoxia (21%O2) or hyperoxia (90% O2) from birth to postnatal day 4 were randomized to receive either hMSCs or vehicle via intratracheal delivery on postnatal day 4. Mouse lungs from normoxia and hyperoxia groups with and without hMSCs therapy were examined at day 14. RNA-sequencing was performed on flow-cytometric CD45+CD11b+CD11c+ sorted lung macrophages. Purified total RNA was used to construct barcoded multiplex-compatible sequencing libraries using: 1) Illumina Stranded mRNA Sample Preparation chemistry (for transcriptomics) and 2) Bio Scientific NEXTFlex Small RNA chemistry (for small RNA). Results Sorted CD45+CD11b+CD11c+ lung macrophages from hyperoxia-exposed neonatal mice showed differentially expressed macrophage genes and miRNA compared to mice exposed to normoxia or hyperoxia+hMSCs. The administration of hMSCs was found to differentially upregulate 421 genes and downregulate 651 genes in CD45+CD11b+CD11c+ lung macrophages from neonatal mice exposed to hyperoxia, compared to normoxia. Integrity pathway analysis (IPA) analysis of macrophage-specific gene pathways revealed the effectiveness of hMSCs in altering macrophage function towards an anti-inflammatory ‘M2’ phenotype. Small-RNA sequencing provided further evidence on the effects of hMSCs, where 1,098 small RNAs transcriptomes were expressed as either significantly up- or down-regulated in response to hMSCs therapy following hyperoxia-induced lung damage. Conclusions Pathway analysis of the predicted mRNA targets of differentially expressed genes provides insight into miRNAs that preferentially target several important pathways. These miRNAs will be functionally relevant for lung macrophages, and will provide a greater understanding of the interaction between macrophage genotype and the associated phenotypes in the setting of inflammation or tissue repair.

Macrophages are traditionally associated with immune functions, but have gained increased recognition for their important roles in tissue regeneration, organogenesis, and mitigation of in ammatory processes [3][4][5]; the phenotype of the macrophages involved in these processes is dependent on the local microenvironment [6]. Lung macrophages are known to play important roles in lung development are have been implicated in the development of lung disease [7], including hyperoxic lung injury in preterm babies [8].
In response to lung injury, immune cell in ltration normally occurs [9]. This includes macrophages with two distinct phenotypes. One is associated with host defence functions, through the production of proin ammatory cytokines and cytotoxic mediators; these are de ned as "classically activated" or M1 [10]. Macrophages also can assume an "alternative" or M2 activation state in which they provide trophic and immunoregulatory functions through tissue remodelling and production of anti-in ammatory mediators [11][12][13]. M2 macrophages can be induced by a variety of stimuli including colony stimulating factor (CSF-1) [14,15]. Lung macrophages are now known to play key regulatory roles in organogenesis and fetal development, as recently reviewed [16,17].
In a recent study [18], we reported the therapeutic effects of human mesenchymal stem cells (hMSCs) on the polarisation of lung macrophages, the development of lung in ammation, and collagen deposition following hyperoxic lung injury in neonatal mice. To date, the therapeutic effectiveness of hMSCs on the gene expression of lung macrophages in response to hyperoxia, and the changes in total and small RNA associated with differing macrophage phenotypes, are yet to be determined. To address this, we have used RNA-Next Generation Sequencing (NGS) to provide transcriptome analysis, the measurement of gene expression, and the identi cation of differential gene regulation in lung macrophages from neonatal lungs exposed to hyperoxia with/without hMSCs therapy [19]. We used NGS diagnostics to explore the early response of lung macrophages to hyperoxia (with and /without hMSCs therapy), and the application of RNA-sequencing that exquisitely de nes the transcriptomes of macrophages in settings of lung in ammation and regeneration.

Experimental animals
The male offspring of C57BL/6J mice obtained from Monash Animal Research Platform (MARP) (ethics approval number MARP2014/092) were used in this study. Pups were born at term, and the day of birth was designated as postnatal day 0 (P0). Throughout the study, the pups were kept in a controlled 12h light and 12h dark lighting cycle with ad libitum access to food and water.
To study the effects of hyperoxic lung injury and its treatment we exposed mouse pups to 90% O 2 from P0 -P4 and used three treatment groups (n=6/group). In the rst group (controls) the pups were raised in room air; the second group were exposed to 90% O 2 ; the third group were exposed to 90% O 2 (as for the second group) and received 2.5x10 5 hMSCs suspended in 10µl PBS delivered via the trachea on P4. Each of the three treatment groups contained six male mice. The injected hMSCs were at passage 4 and purchased from the Tulane Centre for Stem Cell Research and Regenerative Medicine (Tulane University, New Orleans, LA; [20]). At P4, after the period of hyperoxia and hMSCs injection, the pups were kept in room air with their dams until P14 when pups were humanely killed and the lungs collected.
Lung preparation for macrophage isolation At P14, the freshly excised lungs were placed into cold uorescence activated cell sorting (FACS) buffer (PBS containing 0.2% BSA, 0.02% NaN 3 and 5mM EDTA) for a maximum of 30 minutes. The main airways were removed and the lungs were cut into small pieces which were then digested mechanically and enzymatically as previously reported [18].
Lung macrophages identi ed by FACS (CD45 + CD11b + CD11c + ) were sorted via an in ux sorter with a 100 µl nozzle (FlowCore, Melbourne, Australia). The sorted macrophages were washed with PBS, centrifuged for 5 minutes at 2000 rpm at 4°C and then resuspended in carefully aspirated supernatant in 350 µl of RLT lysis buffer with β-mercaptoethanol. Cells were snap-frozen in liquid nitrogen and then stored at -80°C for later RNA extraction as shown in Figure 1.
Total-RNA sequencing of lung macrophages Total RNA was extracted using the RNeasy Kit following the manufacturer's protocol (Qiagen, Hilden, Germany). The quantity and purity of total RNA was evaluated by Nanodrop 1000. From a total of 14µl of extracted RNA, 1µl was used to assess the RNA purity by accepting a Nanodrop ND 1000 reading ratio between 260/280 and 260/230. Total RNA samples were subjected to an initial quality control measurement at Micromon (Monash University-Clayton, Victoria, Australia). They were quantitated using the Invitrogen Qubit and associated chemistry, which incorporates a double-stranded DNA or RNA-speci c uorescent dye (Invitrogen, Carlsbad CA., USA). The integrity of the samples was measured using the Agilent Bioanalyzer 2100 using B.02.08.SI648 (SR3) software, and a micro uidics device, in conjunction with the associated hardware and chemistry (Agilent Technologies, Waldbronn, Germany).
Based on the quality of RNA sample, samples from 3 mice per group (from a group total of n=6) were used to construct Illumina sequencing libraries using the Takara SMART-seq Ultra Low Input RNA kit (version 4), according to the manufacturer's instructions (Clontech Takara SMART-Seq V4 Ultra Low Input RNA Kit for Sequencing User Manual v.01251, California, USA). Libraries were quantitated using a Qubit DNA HS kit, which incorporates a double-stranded DNA-speci c uorescent dye (Invitrogen, Carlsbad CA., USA). The libraries were sized and checked for adapter contamination using the Agilent Bioanalyzer 2100 micro uidics device, in conjunction with Agilent DNA HS kits and chemistry (Agilent Technologies, Waldbronn, Germany). These libraries were sequenced using NextSeq500 -High-Output SBS version 2 (Illumina 15046563 v02), with library concentration 1.8 pico molar (pM) and reading length 1X75b according to the manufacturer's instructions.

Gene expression analysis
A list of differentially expressed genes was generated by comparing the lung macrophage transcriptomes in the samples from the three treatment groups. These genes were further analysed using the Ingenuity Systems Pathway Analysis Knowledge Base (California, USA) to identify underlying cellular gene pathways. Finally, the list of differentially expressed genes generated from analysis of the RNAsequencing data was compared to results from a comparable microarray experiment. One-way Anova was used for statistical analysis.
Macrophage total and small RNA sequencing From each mouse, 1x10 5 CD45 + CD11B + CD11c + live macrophages were sorted from the cells obtained from the whole lung (range 2.7x10 7 -3.1x10 7 ). Genes involved in macrophage maturation [21], activation [22], function [23] and phenotypes [24] have previously been reported. In this project, NGS was used to categorise expressed genes into two sub-categories. Total and small-RNA macrophage gene expression was analysed by assessing the gene integrity pathway [20], which allowed for the analysis of macrophage contribution to hyperoxic lung injury, with and without hMSCs therapy, as reported in our previous study [18]. In the present project, transcriptome analysis of the total RNA allowed us to assess the expression of 11,343 genes. RNA-sequencing was performed by comparing lung macrophages from mice exposed to hyperoxia compared to hyperoxic mice administered hMSCs therapy.

Results
Extracted RNA from macrophages was analysed. with a standard false discovery rate (FDR) of 0.05; most of the genes showed either signi cant up-or down-regulated expression in lungs exposed to 90% O 2 .
Following the hyperoxia, the highest level of gene up-regulation was a (+3.64) Log 2 fold change, found for Cyp4a12b. In contrast, the greatest degree of gene down-regulation was the (-5.21) Log 2 fold change found for Spink1. The administration of hMSCs had an extensive impact on the changes in gene expression level. A total of 421 genes were found to be up-regulated and 651 genes were down-regulated following hMSCs therapy. The changes were expressed as Log 2 fold-changes. In addition, the highest level of gene up-regulation was a (+2.72) Log2 fold change, found for Lpar3, and the lowest level of gene up-regulation was the (+0.30) Log2 fold change, found for Btk. In contrast, a total of 651 genes were down-regulated. The greatest degree of gene down-regulation was the (-4.28) Log2 fold change found for Hspa1b, and the smallest was the (-0.26) Log2 fold change, found for Mapre1.
A total of 139 genes were shown to be up-regulated by hyperoxia and then signi cantly down-regulated following the administration of hMSCs. Of these genes, the highest level of up-regulation was a (+2.74) Log2 fold change, found for Rhov, and the lowest level of gene up-regulation was a (+0.15) Log2 fold change, found for Tcf25, as shown in Figure 2. A list of the 30 most up-regulated genes (among genes upregulated by hyperoxia) following hMSCs administration is shown in Figure 3.

Macrophage total-RNA gene pathways
The use of integrity pathway analysis (IPA) showed the ability of hMSCs to ameliorate the effects of hyperoxia following the analysis of speci c gene expression combinations, that likely contribute to a number of physiological and pathological consequences, as shown in Figure 4. Analysing macrophage-related gene expression and macrophage speci c gene pathways revealed the speci c contribution of macrophages in hyperoxic lung injury, and the effectiveness of hMSCs in altering macrophage function and phenotypes. Using NGS, macrophage total-RNA gene expression was analysed to assess speci c genes characteristic of M1-proin ammatory macrophages, including Cd86, Stat1, Socs3, Slamf1, Tnf, Fcgr1, Il12b, Il6, Il1b, and Il27ra. These genes were up-regulated in lung macrophages isolated from mice exposed to hyperoxia compared to mice receiving hMSCs following hyperoxia, as shown in Figures 5 and 6. Genes related to M2 anti-in ammatory macrophages showed an overall downregulation in macrophages isolated from lungs of hyperoxic mice receiving hMSCs. These genes included Arg1, Stat6, Mrc1, Il27ra, Chil3, and Il12b. Hyperoxia alone led to a signi cant up-regulation in M1associated gene expression together with the down-regulation of M2-associated gene expression. The administration of hMSCs was found to inhibit the expression of M1 macrophage genes and up-regulate M2 associated gene expression. The results are summarised in Figure 7.

Macrophage phenotypes and associated gene pathways
Following the analysis of gene expression, the gene pathways of macrophages were assessed using IPA.
Our total RNA sequencing data showed changes in the expression of macrophage-associated genes. Investigating the contribution of these gene expression changes in macrophages that have associated effects on lung morphology and the development of brosis could clarify the mechanism of the macrophage's response to hyperoxia with and without hMSCs therapy.

Small-RNA sequencing
We used small-RNA sequencing to provide further information on the effects of hyperoxia and the administration of hMSCs on the small RNAs transcriptomes that were evident in the CD45 + CD11b + CD11c + macrophages. The 0.05 FDR was also applied for the small-RNA analysis. A total of 1,098 transcriptomes were expressed as either signi cantly up-or down-regulated following the administration of hMSCs in mice exposed to 90% O 2 . Importantly, these results indicated different types of small RNA transcriptomes expressed that included small nucleolar RNAs (snoRNA and snRNA), ribosomal-RNA (r-RNA), miscellaneous RNA (misc-RNA), microRNA (mi-RNA), large intergenic non-coding RNAs (linc-RNAs) and mitochondrial-RNA (mt-tRNAs).
Up-and down-regulated small-RNA expressiontranscriptomes following hyperoxia A total of 378 small-RNA transcriptomes were up-regulated. The highest-level of transcriptome upregulation was the (+5.04) Log2 fold change found for Gm24144-snRNA and the lowest level of transcriptome up-regulation was the (+0.24) Log2 fold change found for Snord91a-snoRNA. Furthermore, a total of 676 small-RNA transcriptomes were down-regulated. The greatest level of transcriptome downregulation was a (-7.98) Log2 fold change found for Snord88c-snoRNA and the lowest level down transcriptomes down-regulation was (-0.15) Log2 fold change found for Gm24357-snoRNA.
Alterations in transcriptome expression induced by hMSCs in small-RNA transcriptomes that were up and down-regulated by hyperoxia A total of 42 small-RNA transcriptomes were shown to be up-regulated by hyperoxia and then signi cantly inhibited by hMSCs. In these transcriptomes, the highest level of up-regulation was a (+4.71) Log2 fold change, found for A530013C23Rik, and the lowest level of transcriptomes up-regulation was a (+0.56) Log2 fold change, found for Gm22886. In addition, a total of 53 transcriptomes were shown to be down-regulated by hyperoxia and then subsequently signi cantly re-expressed in lungs from mice exposed to hyperoxia and administered hMSCs. The greatest degree of down-regulation in these genes was a (-3.40) Log2 fold change found for Gm24920 and the smallest degree of down-regulation was a (-0.42) Log2 fold change found for Gm13571.

Discussion
Previously we have reported the effects of hMSCs on lung macrophages and their ability to protect neonatal mice from hyperoxia-induced lung injury [18]. Further to this, we have now performed NGS analysis to explore how hMSCs exert their therapeutic effects on macrophage phenotype in established lung injury at the genomic level.
Over the past few decades, various genomic measurement technologies have enabled the gathering of genomic evaluation and the ability to determine whether gene expression is present, absent or modi ed [25]. Technologies such as polymerase chain reaction (PCR) and microarray sequencing [26] use RNA levels to compare gene expression in cells from normal and pathological conditions [27]. However, the degree of novel information offered by NGS is beyond the limits of these traditional RNA and DNA technologies for assessing gene expression.
NGS has been considered by other researchers [28] as a method for studying the genomic basis for hyperoxic lung injury and BPD in preterm babies. We extend these ndings to determine the effect of hMSCs on macrophage phenotype and the resulting improvement of lung structure. In addition, NGS can be used to determine whether these are separate effects, or whether the hMSCs interact with lung macrophages to induce tissue repair. Moreover, NGS can be used to investigate the effects of hMSCs on the genomic expression of other types of immune cells and in related tissue pathologies. hMSCs have remarkable immunomodulatory abilities, and their effects on T cells, B cells, natural killer and dendritic cells have been comprehensively investigated [29,30].
However, the effect of hMSCs on macrophage polarization and the contributions to other cell types that allow for a shift from a pro-in ammatory role to an anti-in ammatory role, in the context of hyperoxic lung injury, remains largely unknown. In this study, by using NGS we have identi ed the effects of hMSCs on shifts in gene expression that correlate with polarization of macrophages towards an M2 phenotype in vivo compared to previous studies that have investigated effects in vitro [31]. We established a successful protocol to isolate macrophage sub-populations using ow cytometry to sort lung macrophages at P14, ten days following the end of the hyperoxia period and the hMSCs injection. We veri ed that at this time-point (P14) there was no evidence of hMSCs remaining in the lung following the initial injection. However, of signi cance, the effects of hMSCs on lung macrophage phenotype were still evident after this extended period.
The extracted total-RNA from the lungs at P14d showed in total, more than ten-thousand genes. It was therefore important to sub-divide these genes into sub-groups, particularly the non-macrophageassociated and macrophage-associated genes. In addition, it was important to exclude the nonsigni cantly changed genes by applying an FDR (false discovery rate) of 0.05. This study supports the rising interest in macrophages that shift their polarisation state and function from pro-in ammatory to anti-in ammatory phenotype, thus contributing to resolving in ammation and improvement of tissue repair after injury [32]. Furthermore, it is known that dysregulation of macrophage polarisation contributes to the development of various in ammatory diseases [33]. The present study provides evidence that hyperoxia upregulated the expression of several genes including Cyp4a12b Tnc, Fbln2, Cyp4a12a, Cd209e, Alas2, Rhov, Spic among others.
Following the identi cation of a wide range of genes of interest, the major new ndings of this study extend to the analysis of transcriptomic changes in hyperoxic lung injury. These results showed that there are key regulators that stimulate expression of lung injury genes, in addition to being involved in injury and in ammation in other organs. This suggests that lung damage might not only develop directly from hyperoxia per se, but also indirectly via altered expression of genes that are associated with tissue damage. Other studies have shown a positive correlation between hyperoxia and abnormal growth of cardiac endothelial cells [34]. Furthermore, hyperoxia is also involved in the development of hepatic pathology, including liver in ammation. Moreover, hyperoxia causes a signi cant aggravation of hepatic tissue damage [35].
Finally, we used the pathway analysis tool to identify the functions and the relations of lung macrophageassociated gene expression. Based on a well-established list of genes and transcriptional patterns that categorise macrophage phenotypes, this study provides evidence of the ability of hMSCs to alter macrophage gene expression by inhibiting the up-regulation of M1-macrophage-associated genes and up-regulating the M2-macrophage-associated genes. Previous studies have shown the transcriptional pattern of lung macrophages and the pattern of macrophage polarization-related genes that are classi ed into M1-and M2-related genes [36][37][38]. Furthermore, our results are compatible with other ndings that have established the ability of hMSCs to alter macrophage polarisation towards an M2 phenotype [39][40][41].
These results showed that the use of IPA facilitates the identi cation of the genes involved with an M2 phenotype, including IL6, STAT6, MRC1 and ARG1, known for their effective anti-in ammatory and regenerative abilities, as well as other genes that are associated with protection against lung injury and brosis in both mice and humans [42].
Together, our studies highlight the ability of hMSCs to target many factors involved in hyperoxia-induced lung injury. Their ability to reduce brosis, repair damaged tissue and improve regenerative environments in the lung make hMSCs a prime candidate for a cell-based therapy to reduce hyperoxia-induced lung injury. Furthermore, the immunosuppressive properties of hMSCs, which can reduce the in ammatory environment, make them an ideal therapy to reduce the symptoms of BPD.

Conclusions
In summary, we report that hyperoxia in neonatal mice increases gene expression of M1-macrophages with an associated decrease in M2 macrophage genes. hMSC-derived factors were shown to impair these changes in gene expression and promote the polarisation towards an alternatively activated M2 phenotype, as shown by their gene expression pro le and gene pathways. In preterm babies, reducing and controlling in ammation is vital for their survival. Therefore, the use of hMSCs could serve as a useful tool for the development of novel therapeutic options for preterm babies with BPD and other respiratory illnesses.

Consent for publication
Not applicable Availability of data and materials The datasets generated during and/or analyzed during the current study are available from the corresponding author on reasonable request.

Competing interests
The authors declare that they have no competing interests.  Figure 1 Diagram showing a summary of the experimental protocol. Neonatal mice were exposed to either 90% O2

Figures
(hyperoxia) or 21% O2 (normoxia, room air) from P0d to P4d. At P4d the hyperoxia was discontinued and a sub-group from the hyperoxia group received 2.5x105 hMSCs injected intratracheally. At P14d all mice were killed and their lungs removed (n=6/group). Lung macrophages were isolated and macrophage RNA was extracted for analysis by next generation sequencing (NGS). Finally, gene expression and pathways were analysed. Abbreviations: hMSCs, human mesenchymal stem cells; P4d, postnatal day 4; and RBC, red blood cell. Up-regulated and down-regulated genes in lung macrophages following hyperoxia in comparison to hyperoxia-treated mice following hMSCs therapy. (A) Heatmap showing gene expression from total-RNA extracted from lung macrophages of mice exposed to hyperoxia, compared to mice with hyperoxiainduced lung injury administered hMSCs. (B) Bar chart showing genes that were increased in response to hyperoxia and altered gene expression in response to hMSCs -total of 139. The greatest differentially expressed gene correlated to Rhov and the least expressed was identi ed as Tcf25. (C) Bar chart showing genes that were down-regulated in expression in lung macrophages from hyperoxic mice that were inhibited via delivery of hMSCs -total of 99. The greatest value for downregulated expression was found in Rab39 and the smallest value was found in Rgs19. Both of the heatmap and the gene expression were measured on log2 fold of change.

Figure 3
Table listing the 30 most up-regulated genes in lung macrophages at P14d following exposure to hyperoxia that were inhibited by hMSCs therapy. These genes were upregulated by hyperoxia alone (90% O2) administered from P0-P4. The table columns denote (H) hyperoxia; (H+hMSCs) mice exposed to hyperoxia and receiving hMSCs, and (FDR) false discovery rate. Data for (Log2 fold-change) are the average of one sample from each of 3 mice in each treatment group. Gene expression (Log2 fold change values) in lung macrophages from hyperoxia-exposed mice with/without hMSCs therapy. Gene expression in lung macrophages from mice exposed to either normoxia, hyperoxia or following the administration of hMSCs (n=3/group). The gure shows the Log2fold change in gene expression up-regulated by hyperoxia and showing a signi cant reduction towards baseline following hMSCs delivery (false discovery rate FDR=0.05, *P-value <0.05) for each expression.   Table shows showing the Log2 fold change value of macrophage genes (M1 and M2 phenotypes). Gene expression in lung macrophages from mice exposed to either normoxia, hyperoxia or following the administration of hMSCs (n=3/group). The genes that are characteristic of M1 and M2 macrophages were either up-regulated or down-regulated in response to hyperoxia. However, the administration of hMSCs to hyperoxic mice was able to alter the expression towards baseline levels. The  M1 macrophage-associated gene expression (Log2 fold change values) in lung macrophages from hyperoxia-exposed mice with/without hMSCs therapy. The graphs show the Log2-fold change in gene expression that denote M1 macrophages up-regulated by hyperoxia and showed a signi cant reduction towards baseline following hMSCs delivery (false discovery rate FDR=0.05, *P-value <0.05) for each expression. Abbreviations: Cd86, Cluster of Differentiation 86; Socs3, suppressor of cytokine signalling 3; Tnf, tumour necrosis factor and Il12b, Interleukin 12B; and hMSCs, human mesenchymal stem cells.