Cytochrome P450 26A1 regulates the clusters and killing activity of NK cells during the peri‐implantation period

Abstract Cytochrome P450 26A1 (CYP26A1) plays a vital role in early pregnancy in mice. Our previous studies have found that CYP26A1 affects embryo implantation by modulating natural killer (NK) cells, and that there is a novel population of CYP26A1+ NK cells in the uteri of pregnant mice. The aim of this study was to investigate the effects of CYP26A1 on the subsets and killing activity of NK cells. Through single‐cell RNA sequencing (scRNA‐seq), we identified four NK cell subsets in the uterus, namely, conventional NK (cNK), tissue‐resident NK (trNK) 1 and 2, and proliferating trNK (trNKp). The two most variable subpopulations after uterine knockdown of CYP26A1 were trNKp and trNK2 cells. CYP26A1 knockdown significantly downregulated the expression of the NK cell function‐related genes Cd44, Cd160, Vegfc, and Slamf6 in trNK2 cells, and Klra17 and Ogn in trNKp cells. Both RNA‐seq and cytotoxicity assays confirmed that CYP26A1+ NK cells had low cytotoxicity. These results indicate that CYP26A1 may affect the immune microenvironment at the maternal‐foetal interface by regulating the activity of NK cells.


| INTRODUC TI ON
Pregnancy is a unique biological event in which a mother peacefully coexists with a semi-allogeneic foetus. During the first trimester of pregnancy, leukocytes account for ~30%-40% of cells in the human decidua. Natural killer (NK) cells make up ~70% of decidual leukocytes. 1,2 Uterine NK cells in mice are defined as CD45 + Lin − NK1.1 + NKp46 + during pregnancy. Two different functional NK cell populations can be identified by their ability to bind to the lectin Dolichos biflorus agglutinin (DBA) in mice. 3 DBA + NK cells lack the expression of CD49b but overwhelmingly express angiogenic factors, whereas Ifng expression prevails among the DBA − NK cells which are CD49b positive. 4,5 Tian et al. 6 define CD49a as a specific cell surface marker of NK cells residing in tissue. Based on this, murine NK cells are classified into CD49a − conventional NK (cNK) cells and CD49a + tissue-resident NK (trNK) cells. 7,8 trNK cells are dominant during early gestation, most abundant on gestation day (GD)5.5, and decrease in number as pregnancy progresses. After placenta formation, cNK cells expand rapidly and dominate during late pregnancy. 9 NK cells, the major leukocyte population at the maternal-foetal interface, play a vital role in embryo implantation, trophoblast invasion and spiral artery remodeling during pregnancy. 10,11 NK cell exhaustion in pregnant mice leads to a reduction in the number of implantation sites, barren vascular remodelling and an increase in the embryo resorption rate. 12,13 In humans, NK cells infiltrate the decidua and locate nearby extravillous trophoblast cells (EVTs) during early pregnancy. Certain combinations of human leukocyte antigen C expressed by EVTs and killer immunoglobulin-like receptors on decidual NK (dNK) cells contribute to increased risk of preeclampsia. This is due to excessive inhibition of dNK cells, which results in poor trophoblast invasion. 14 Conversely, unduly active uterine NK cells are closely related to pregnancy failure. 15,16 Cytochrome P450 26A1 (CYP26A1), a member of the cytochrome P450 superfamily, is a monooxygenase that catalyses the metabolism of all-trans-retinoic acid (at-RA). 17 Cyp26a1-null mice die mid-pregnancy or during terminal pregnancy and exhibit numerous crucial morphogenetic defects, such as aberrant hindbrain patterning and vertebral identity. 18 Our previous work has found that CYP26A1 shows a peculiar temporal and spatial expression pattern in mice and rats during the peri-implantation period. 19,20 In Cyp26a1-MO-treated and pCR3.1-Cyp26a1 plasmid-immunized mice, the number of implantation sites significantly decreases and the proportion of NK cells, dendritic cells (DCs) and macrophages changes dramatically. [21][22][23] dNK1 cells that highly express CYP26A1 have been identified during early gestation in humans. 24 Previous studies in our laboratory have also found a novel population of CYP26A1 + NK cells in the uterus. 25 However, little is known about the underlying mechanisms involved, especially the effect of CYP26A1 on NK cell immune activity.
In this study, single-cell RNA sequencing (scRNA-seq) analysis revealed four major NK cell subsets in the uterus. CYP26A1 knockdown had no effect on the clusters of NK cells but it affected the proportion of NK cell subsets and downregulated the expression of genes associated with immune activity and cytokines in NK cells.
Further experiments indicated that CYP26A1 + NK cells had low killing activity. Hence, CYP26A1 may influence the immune microenvironment at the maternal-foetal interface by modulating the activity of NK cells.

| Mice
BALB/c mice, aged 8-10 weeks, were purchased from Vital River Laboratory Animal Technology Co. Ltd. (Beijing). The mice were housed in the animal care facility of the Institute of Zoology, Chinese Academy of Sciences, according to the institutional guidelines for the care and use of laboratory animals. Male mice cohabited with females at a ratio of 1:1. The day when the vaginal plug was detected was recorded as the first day of pregnancy (GD1).

| Preparation of single-cell suspension and flow cytometry analysis
The uterus and spleen were cut apart and minced into small fragments. Splenic tissue fragments were placed in PBS containing 2% FBS. Uterine fragments were placed in 1640 medium containing 200 U/ml hyaluronidase (H3506, Sigma-Aldrich), 1 mg/ml collagenase type IV (C5138, Sigma-Aldrich) and 8% FBS and then incubated at 37°C for 30 min, as previously described but with minor changes. 26 After digestion, uterus cells were centrifuged and incu-

| Western blotting
Mouse uterine tissue was ground into powder in liquid nitrogen and added to RIPA Lysis Buffer (CW2334S, Cwbio) containing 1 mM PMSF (78830, Sigma). Total protein was extracted following the instructions of the RIPA Lysis Buffer kit. The Bicinchoninic Acid Protein Assay Kit (23227, Pierce) was used to detect the protein concentration. Proteins were separated using 10% SDS-PAGE and then transferred from the gel onto a nitrocellulose membrane (66485, Pall). The membrane was blocked with 5% skimmed milk and then incubated overnight with primary antibodies at 4°C. After washing thoroughly with TBST solution, the membranes were incubated with HRP-coupled secondary antibodies at room temperature for 1 h and then visualized using a chemiluminescence imaging system (MiniChemi 610, Sagecreation). The data were analysed using ImageJ software. The primary and secondary antibodies used for Western blotting included anti-CYP26A1 (ab151968, 1:1000, Abcam), anti-GAPDH (2118, 1:1000, Cell Signaling Technology) and goat antirabbit IgG (H + L) HRP (31460, 1:10,000, Thermo Fisher).

| scRNA-seq and analysis
Viable CD45 + CD3 − CD122 + cells were sorted from the uteri of five Cell Ranger software was used to process raw sequencing data, and the mm10 mouse reference genome assembly was used as the reference index. Cells were removed when one of the following conditions was met: (i) the number of detected genes was more than 8000 or less than 300; (ii) the number of UMIs was less than 500; and (iii) the mitochondrial gene expression level was higher than 20%. In addition, genes expressed in more than three cells were kept.
Raw count matrices, obtained by reading the output of the Cell Ranger pipeline, were converted into Seurat objects, which were merged together into a single Seurat object for downstream analysis according to Seurat (version 4.0.6) tutorials. 28 These analyses, including normalization, variance stabilization, sample integration using shared hypervariable genes, cell clustering based on top principal components, differential expression analysis and visualization, were carried out following the standard process of Seurat. To characterize the functional properties of four NK cell subsets, KEGG analysis of differentially expressed genes (Absolute log 2 Fold change > 0.5 and adjusted p < 0.05) among these NK cell subsets was performed using the R package clusterProfiler (version 4.2.0). 29 clusterProfiler was also used for GO analysis of differentially expressed genes (adjusted p < 0.05; |log 2 Fold change| > 0.5) among Cyp26a1-MO treated mice and controls.

| Immunocytofluorescence
NK cells were isolated from GD5 mice by flow cytometry. After centrifugal re-suspension, cells were incubated with DAPI (2 μg/ml) at room temperature for 15 min. After washing, the cell suspensions were dropped onto a glass slide and covered with a cover slip for imaging. Images were captured using a Zeiss LSM 780 confocal microscope and analysed with ZEN software.
They were then co-cultured with effector NK cells at an effectorto-target cell (E:T) ratio of 2:1 or 0:1 for 4 h at 37°C under 5% CO 2 .
Fluorescence-activated cell sorting (FACS) analysis was performed immediately with a BD LSR Fortessa Cell Analyzer (BD Biosciences).
The Incucyte S3 Live-Cell Analysis Instrument (Sartorius AG) was used to monitor the killing activity of NK cells. Propidium iodide (C1052, Beyotime) was used to mark dead cells.

| Statistical analysis
Data were analysed with GraphPad Prism (version 9.0.0) software.
The results were expressed as means ± SEM. The paired or unpaired two-tailed t-test was used to evaluate the differences. Statistically significant differences were defined as p values < 0.05.

| Analysis of NK cell clustering via scRNA-seq
To and then performed scRNA-seq using the 10x platform ( Figure 1C).
After filtering out the low-quality cells, we finally acquired the tran- among the three trNK cell subsets and were characterized by Aff3 ( Figure 1E). CYP26A1 knockdown had no effect on the clustering of uterine NK cells.

| CYP26A1 affected NK cell proportion and activity
As mentioned above, we defined four main NK cell subsets: cNK, trNKp, trNK1 and trNK2. To determine the functional characteristics of these subsets, markers were found for each cluster using Seurat; the expression levels of the top 10 marker genes in each subset are presented in Figure 2A. KEGG pathway analysis was performed using the differentially expressed marker genes (adjusted p < 0.05; |log2Fold change| > 0.5) for each subset ( Figure 2B). As shown, cNK Klra17 and the growth-promoting factor Ogn after CYP26A1 knockdown ( Figure 2E). The activation of immune response pathway was enriched in trNKp cells ( Figure 2F).

| CYP26A1 + NK cells underwent dynamic changes
To verify the expression of CYP26A1 in NK cells, we employed flow cytometry to detect the proportions of CYP26A1 + NK cells isolated from mouse spleens and uteri on GD6. Among mouse uterine cells, approximately 9.5% of CD45 + CD3 − cells were CD122 + CYP26A1 + cells, whereas only 0.6% spleen cells were CD122 + CYP26A1 + cells ( Figure 3A,B). Flow cytometry analysis indicated that CYP26A1 + NK cells had a higher median fluorescence intensity (MFI) value for CYP26A1 than CYP26A1 − cells ( Figure 3C). Live-cell immunofluorescence also showed that CYP26A1 was only present on CYP26A1 + NK cells and was primarily located on the cytoplasmic membrane ( Figure 3D). These data demonstrated that there was a specific CYP26A1 + NK cell subset in the uteri of pregnant mice.
Before investigating the impact of CYP26A1 + NK cells on early pregnancy, we detected the changes in the trend of their proportion during the peri-implantation period. DBA, which has a high affinity for glycoconjugates containing N-acetyld-galactosamine in the terminal position, is known as the specific marker for mouse uterine NK cells. 3 Flow cytometry was used to analyse the proportions of CYP26A1 + and DBA + cells among NK cells marked as CD45 + CD3 − CD122 + in mice from GD4 to GD7 ( Figure 3E). The proportion of CYP26A1 + cells was highest (34.10%) on GD4, decreased dramatically on GD5 (p < 0.0001), and was maintained at a relatively stable level (~20%) from GD5 to GD7 ( Figure 3F). The proportion of DBA + cells peaked on GD4 (6.90%), decreased significantly from GD4 to GD5 (p = 0.0011), and showed a further significant reduction from GD6 to GD7 (p = 0.0138) ( Figure 3G).

| Killing activity of CYP26A1 + NK cells
To investigate the function of CYP26A1 + NK cells, we performed a flow cytometry-based cytotoxicity assay to evaluate the killing activity of freshly isolated splenic and uterine NK cells on GD6.
Briefly, YAC-1 targets were labelled with CFSE and then co-cultured with effector NK cells at an E:T ratio of 2:1 for 4 h ( Figure 4A). In the assay, non-viable apoptotic and dead cells were marked as Annexin V + NIR + and NIR + cells respectively ( Figure 4B). Uterine CYP26A1 + NK cells showed the lowest percentage of non-viable apoptotic and dead cells within YAC-1 cells, followed by uterine and splenic CYP26A1 − cells ( Figure 4C). Further experiments were conducted to assess the lytic capacity of uterine CD49b + , CD49b − , CYP26A1 − and CYP26A1 + NK cells. The results indicated that CYP26A1 + NK cells had lower killing capacity than the other three uterine NK cell subsets ( Figure 4D). We also monitored cell viability in real-time using an Incucyte S3 live-cell analysis system. YAC-1 cells co-cultured with uterine CYP26A1 + NK cells showed a higher viability than CYP26A1 − NK cells from the spleens and uteri of pregnant mice ( Figure 4E).
These results demonstrated that uterine CYP26A1 + NK cells had low killing activity during early pregnancy.

| CYP26A1 + NK cells with specifical transcriptional profile
To further investigate functional differences between uterine CYP26A1 − and CYP26A1 + NK cells in mice during the peri-implantation period, we conducted comprehensive transcriptome-wide screening to evaluate their transcriptional expression profiles. Transcriptome analysis revealed 4612 differentially expressed genes (adjusted p < 0.05; |Fold change| > 2). Among them, 1,669 genes displayed lower expression, whereas 2,943 genes were upregulated in CYP26A1 + NK cells ( Figure 5A). We then used these differentially expressed genes to perform GSEA. The 10 KEGG pathways significantly enriched in GSEA are shown in Figure 5B.
We also performed gene overlap relationship analysis on the top 5 pathways enriched by KEGG ( Figure 5C). Interestingly, the natural killer cell-mediated cytotoxicity pathway was enriched with an adjusted p value of 0.001 ( Figure 5B,C). The pathway was dramatically downregulated in CYP26A1 + NK cells ( Figure 5D). Compared with CYP26A1 − NK cells, CYP26A1 + NK cells displayed lower expression levels of pivotal genes involved in this pathway, such as Prf1, Gzmb and Fasl ( Figure 5E). Activating and inhibitory receptors also modulated NK cell activity. We found that the majority of activating receptors, such as Cd226, Klrk1 and Ncr1 had lower expression levels in CYP26A1 + NK cells ( Figure 5F). The levels of inhibitory receptors such as Siglece, Fcgr2b and Klra2 were higher in CYP26A1 + NK cells ( Figure 5G).

| RNA-seq data validation
Using GSEA, we found that the natural killer cell-mediated cytotox-

| DISCUSS ION
CYP26A1, an RA-metabolizing enzyme, has already been verified to play a prominent role in embryo implantation. Both inhibition and blockade of CYP26A1 function results in a significant reduction in the number of embryo implantation sites in mice. 19 Further studies have indicated that the at-RA concentration in the uterus has no marked difference after CYP26A1 knockdown (data not published) and that intraperitoneal administration of supraphysiological doses of at-RA has no influence on embryo implantation in mice. 23 These results indicate that CYP26A1 may regulate embryonic implantation via a non-RA pathway. Recently, we found that CYP26A1 regulates the differentiation of DCs (through CD86 and ID2), polarization of uterine macrophages, and proportion of NK cells during the periimplantation period in mice. [21][22][23] We hence conclude that CYP26A1 affects embryo implantation through immune cells at the maternalfoetal interface.
To the best of our knowledge, this study is the first to provide a single-cell transcriptomics atlas of NK cells at the maternal-foetal interface in Cyp26a1-MO knockdown mice. We defined four major  However, the exact mechanism needs to be clarified through further experiments.
Our laboratory has found that there is a population of CYP26A1 + NK cells at the maternal-foetal interface. 25 Unfortunately, we detected only five NK cells expressing Cyp26a1 via scRNA-seq (data not shown). The most likely explanation is that the latter had limited sequencing depth. Flow cytometry and immunofluorescence assays confirmed that the CYP26A1 + NK cell subset specifically existed in the uteri of pregnant mice, which dovetails nicely with the results for humans. 24,35 This NK cell subset exhibited dynamic changes during the peri-implantation period. Cytotoxicity assays demonstrated that CYP26A1 + NK cells had low killing activity. Given that NK cell activ- In conclusion, our data indicated that CYP26A1 knockdown had no effect on the clusters of NK cells in the uterus, but it altered the proportion of NK cell subsets and significantly downregulated the expression levels of cytokines and immunologic activity-related genes in NK cells. CYP26A1 + NK cells exhibited low killing activity. CYP26A1 may affect the immune microenvironment at the maternal-foetal interface by regulating the activity of NK cells.

ACK N OWLED G EM ENT
This work was supported by grants from the National Natural Science Foundation of China (no. 31571552).

CO N FLI C T O F I NTE R E S T
The authors confirm that there are no conflicts of interest.

DATA AVA I L A B I L I T Y S TAT E M E N T
We confirm that the data in our paper can be used.