Alternative splicing and differential gene expression during changes in endometrial receptivity in patients with recurrent implantation failure.

The patient recruitment, sample collection, and transcriptomic data analysis of healthy controls have been previously described 6 . No statistically significant differences were observed in age, body mass index, previous implantation failure, menstrual cycle length, or endometrial thickness among the three groups of endometrial samples collected from RIF patients undergoing HRT (Table 1 ). The average raw RNA-seq data generated from endometrial tissue samples of RIF patients was approximately 9.70 GB. After quality control and removal of ribosomal RNA (see Materials and Methods), the average clean data size was 8.89 GB. All clean reads were mapped to the human reference genome (GRCh38), and each sample achieved a mapping rate exceeding 80%. To obtain the landscape of alternative splicing variation across human endometrial receptive phases, we performed differential expression and alternative splicing analyses between the non-receptive (PR or PS) and receptive (R) phases in both healthy controls and RIF patients. In total, 1,240,563 ASEs were identified between PR and R stages in healthy controls (HC: PR vs. R), 981,422 ASEs between PS and R stages in healthy controls (HC: PS vs. R), 985,440 ASEs between PR and R stages in RIF patients (RIF: PR vs. R), and 1,010,588 ASEs between PS and R stages in RIF patients (RIF: PS vs. R), with percentage spliced-in (PSI) scores calculated (Table 2 ). Table 2 Statistical results of alternative splicing and differential expression analyses in four different subgroups. Group EventType ALL Events Differential ASEs Differential ASEs (Up) Differential ASEs (Down) ALL Events (Total) Differential ASEs (Total) ASEGs (Total) DEGs (Up) DEGs (Down) DEGs (Total) HC PR vs. R SE 491,214 520 229 291 1,240,563 1,048 505 600 650 1,250 A5SS 168,616 88 48 40 A3SS 233,885 114 75 39 MXE 285,928 239 112 127 RI 60,920 87 59 28 PS vs. R SE 403,765 276 150 126 981,422 544 342 258 232 490 A5SS 138,542 35 13 22 A3SS 190,273 52 26 26 MXE 195,983 127 43 84 RI 52,859 54 24 30 RIF PR vs. R SE 406,658 911 483 428 985,440 1,744 792 1,078 1,988 3,066 A5SS 136,024 126 90 36 A3SS 195,840 155 118 37 MXE 192,212 376 165 211 RI 54,706 176 143 33 PS vs. R SE 418,804 431 239 192 1,010,588 804 437 394 652 1,046 A5SS 137,771 53 32 21 A3SS 200,145 78 38 40 MXE 198,746 161 68 93 RI 55,122 81 65 16 HC: healthy control; RIF: recurrent implantation failure; PR: pre-receptive; R: receptive; PS: post-receptive; SE: skipped exons; A5SS: alternative 5’ splice sites; A3SS: alternative 3’ splice sites; MXE: mutually exclusive exons;RI: retained introns; ASEs: alternative splicing events; ASEGs: genes with differentially alternative splicing event; DEGs: differentially expressed genes. Statistical results of alternative splicing and differential expression analyses in four different subgroups. HC: healthy control; RIF: recurrent implantation failure; PR: pre-receptive; R: receptive; PS: post-receptive; SE: skipped exons; A5SS: alternative 5’ splice sites; A3SS: alternative 3’ splice sites; MXE: mutually exclusive exons;RI: retained introns; ASEs: alternative splicing events; ASEGs: genes with differentially alternative splicing event; DEGs: differentially expressed genes. After applying stringent criteria (FDR  0.1), 1,048 differential ASEs were identified in the HC: PR vs. R comparison, 544 in HC: PS vs. R, 1,744 in RIF: PR vs. R, and 804 in RIF: PS vs. R (Table 2 ; Figure S1 A,B). Consistent across all groups, skipped exons (SE) and mutually exclusive exons (MXE) were the most common splicing events. Further analysis of the PR and R phases in healthy controls revealed 1,048 differential ASEs, including SE (n = 520, 49.62%), MXE (n = 239, 22.80%), A5SS (n = 88, 8.40%), A3SS (n = 114, 10.88%), and RI (n = 87, 8.30%). In the RIF: PR vs. R group, 1,744 ASEs were identified, comprising 52.24% SE, 21.56% MXE, 7.22% A5SS, 8.89% A3SS, and 10.09% RI (Fig.  2 A). Fig. 2 Alternative splicing and immune cell analyses between pre-receptive and receptive phases in RIF patients. (A) Analysis of differential ASEs between endometrial PR vs. R phases of HC and RIF (green) patients (orange). The percentage of each type of ASE within each comparison is represented by a pie chart (right) (SE, skipped exon; MXE, mutually exclusive exons; A3SS, alternative 3’ splice site; A5SS, alternative 5’ splice site; RI, retained intron). (B) Venn diagram showing the genes with alternative splicing events in HC: PR vs. R (blue) and RIF: PR vs. R (orange). (C) The proportions of differential ASEs for each type in in HC: PR vs. R (blue) and RIF: PR vs. R (orange). (D) The frequency of possible changes in coding potential shown for SE, MXE, and RI, respectively. (E) Scatter plot displayed the GO analysis of genes with different alternative splicing events. BP: Biological process, CC: cellular component, MF: molecular function. (F) Comparing immune cell distribution between endometrial PR vs. R phase of RIF patients was shown in the boxplot (ns: p > 0.05, * p < = 0.05, ** p < = 0.01, *** p < = 0.001, **** p < = 0.0001). (G) Heatmap of the distribution of immune cells interacting with different ASEs in the endometrial R phases of RIF patients. Alternative splicing and immune cell analyses between pre-receptive and receptive phases in RIF patients. (A) Analysis of differential ASEs between endometrial PR vs. R phases of HC and RIF (green) patients (orange). The percentage of each type of ASE within each comparison is represented by a pie chart (right) (SE, skipped exon; MXE, mutually exclusive exons; A3SS, alternative 3’ splice site; A5SS, alternative 5’ splice site; RI, retained intron). (B) Venn diagram showing the genes with alternative splicing events in HC: PR vs. R (blue) and RIF: PR vs. R (orange). (C) The proportions of differential ASEs for each type in in HC: PR vs. R (blue) and RIF: PR vs. R (orange). (D) The frequency of possible changes in coding potential shown for SE, MXE, and RI, respectively. (E) Scatter plot displayed the GO analysis of genes with different alternative splicing events. BP: Biological process, CC: cellular component, MF: molecular function. (F) Comparing immune cell distribution between endometrial PR vs. R phase of RIF patients was shown in the boxplot (ns: p > 0.05, * p < = 0.05, ** p < = 0.01, *** p < = 0.001, **** p < = 0.0001). (G) Heatmap of the distribution of immune cells interacting with different ASEs in the endometrial R phases of RIF patients. Among the 505 ASEGs detected in the HC: PR vs. R group, only 26 (5.15%) were also DEGs, whereas 479 (94.85%) were not differentially expressed at the gene level. Similarly, in the RIF: PR vs. R group, only 93 (11.74%) of the 792 ASEGs overlapped with DEGs (Figure S2 A). A total of 178 ASEGs were shared between the HC: PR vs. R and RIF: PR vs. R comparisons, of which 11 genes ( CKB, FN1, GRAMD1C, KIF12, NNMT, PDGFA, RABGAP1L, RIMKLB, SDCBP2, SYNE2, TYMP ) were also differentially expressed (Fig.  2 B; Figure S2 A). Among these, 31 ASEGs have been previously reported as associated with hormonal regulation, endometrial stromal cell differentiation, trophoblast development, and implantation (Table 3 ). Five of these genes ( EPB41L2, OFD1, POLD4, POSTN, and SH3YL1 ) have been identified as predictive markers of the endometrial WOI in RIF patients 7 , 12 . Six other differentially expressed ASEGs ( CKB, GRAMD1C, KIF12, NNMT, PDGFA, and SYNE2 ) were reported to be related to endometrial receptivity 7 , 12 , 46 , 47 . Table 3 List of 43 overlapping genes with differential alternative splicing events associated with embryo implantation. Number Gene symbol Gene name Relevant function or diseases description PMID 1 ACTN1 a,b a,b Actinin Alpha 1 Regulation by leukaemia inhibitory factor in uterine luminal epithelial cells 25,031,358 2 AIMP1 a,b a,b Aminoacyl TRNA Synthetase Complex Interacting Multifunctional Protein 1 Receptivity associated gene 23,555,582 3 ANXA1 a,b a,b Annexin A1 Regulation of steroid hormone secretion; Maintenance of the uterine microenvironment during implantation 32,403,233;29,115,663;22,819,759 4 ANXA2 b b Annexin A2 Embryo adhesiveness 22,645,245;33,010,173 5 APOL2 b,d b,d Apolipoprotein L2 Decidualization 20,008,415 6 C3 b b Complement C3 Genes for endometrial receptivity prediction 20,619,403 7 CD44 a,b a,b CD44 molecule (Indian blood group) Dynamically expressed across the menstrual cycle; Endometrial stromal cell proliferation and decidualization; Embryo adhesion 8,560,955;16,932,025;36,527,033 8 CKB a,c a,c Creatine Kinase B Genes for endometrial receptivity prediction 20,619,403 9 COX6C b b Cytochrome C Oxidase Subunit 6C Early embryo invasion 34,643,467 10 CREM b b CAMP Responsive Element Modulator Decidualization 21,159,852 11 CTSB a a Cathepsin B Embryonic development and Endometrial metamorphosis 9,310,336 12 DCN a a Decorin Receptivity associated gene 23,555,582 13 DLX6-AS1 a a DLX6 Antisense RNA 1 Dynamically expressed across the menstrual cycle; Regulation of trophoblast proliferation, migration and invasion in patients with pre-eclampsia 28,395,321;30,055,134 14 ELP3 b,d b,d Elongator Acetyltransferase Complex Subunit 3 Embryonic development 27,476,491 15 EPB41L2 a a Erythrocyte Membrane Protein Band 4.1 Like 2 Genes for endometrial receptivity prediction 33,910,562 16 FAP a a Fibroblast Activation Protein Alpha Receptivity associated gene 23,555,582 17 FN1 a,c a,c Fibronectin 1 FN1 isoforms and receptors associated with preimplantation embryo development 19,126,199 18 GABRP a a Gamma-Aminobutyric Acid Type A Receptor Subunit Pi Human early pregnancy trophoblast markers 32,359,161 19 GCH1 b b Embryo lethality 25,557,619 20 GRAMD1C a,c a,c GRAM domain containing 1C Receptivity associated gene 23,555,582 21 GSN a a Gelsolin Epithelial remodeling and embryo adhesion 29,763,784 22 HMGN1 a a High Mobility Group Nucleosome Binding Domain 1 Decidualization of uterine stromal cells 26,566,865 23 HNRNPH1 a a Heterogeneous Nuclear Ribonucleoprotein H1 Regulation of alternative splicing in germ cells; male infertility 35,739,118 24 IL1R1 a a Interleukin 1 Receptor Type 1 Receptivity associated gene 23,555,582 25 KIF12 a,c a,c Kinesin Family Member 12 Genes for endometrial receptivity prediction 33,910,562 26 NNMT a,c a,c Nicotinamide N-Methyltransferase Genes for endometrial receptivity prediction 33,910,562;20,619,403 27 NR1H3 a a Nuclear Receptor Subfamily 1 Group H Member 3 Human trophoblast invasion 15,242,983;18,276,933 28 NUCB2 b b Nucleobindin 2 Ovarian steroidogenesis and uterine function local regulator 30,981,497 29 OFD1 a a OFD1 Centriole And Centriolar Satellite Protein Genes for endometrial receptivity prediction 33,910,562;20,619,403 30 PAX8 a a Paired Box 8 Receptivity associated gene 23,555,582 31 PDGFA a,c a,c Platelet Derived Growth Factor Subunit A Receptivity associated gene 25,429,785 32 PLXNB2 a a Plexin B2 Integrity of endometrial epithelium 25,237,006 33 POLD4 a a DNA Polymerase Delta 4, Accessory Subunit Genes for endometrial receptivity prediction 20,619,403 34 POSTN a a Periostin Genes for endometrial receptivity prediction 20,619,403 35 SEC61A1 b b SEC61 Translocon Subunit Alpha 1 Decidualization 32,386,616 36 SECISBP2L b b SECIS Binding Protein 2 Like Embryonic development 35,210,313 37 SGK1 a a Serum/Glucocorticoid Regulated Kinase 1 Receptivity associated gene;maintenance of pregnancy 22,001,908;27,871,060 38 SH3YL1 a a SH3 And SYLF Domain Containing 1 Genes for endometrial receptivity prediction 33,910,562 39 SLC3A2 b b Solute Carrier Family 3 Member 2 Trophoblast differentiation 35,273,963 40 SYNE2 a,c a,c Spectrin Repeat Containing Nuclear Envelope Protein 2 Genes for endometrial receptivity prediction 20,619,403 41 TJP1 a a Tight Junction Protein 1 Human trophoblast proliferation and invasion 35,687,903 42 TSC2 b b TSC Complex Subunit 2 Follicular depletion and low fertility in female mice 22,128,018 43 UCA1 a a Urothelial Cancer Associated 1 Endometrial stromal cell autophagy and apoptosis 33,680,939 a: Genes with differentially alternative splicing events that intersect in groups HC:PR vs. R and RIF:PR vs. R. b: Genes with differentially alternative splicing events that intersect in groups HC:PS vs. R and RIF:PS vs. R. c: DEGs intersect in groups HC:PR vs. R and RIF:PR vs. R.d: DEGs intersect in groups HC:PS vs. R and RIF:PS vs. R. List of 43 overlapping genes with differential alternative splicing events associated with embryo implantation. a: Genes with differentially alternative splicing events that intersect in groups HC:PR vs. R and RIF:PR vs. R. b: Genes with differentially alternative splicing events that intersect in groups HC:PS vs. R and RIF:PS vs. R. c: DEGs intersect in groups HC:PR vs. R and RIF:PR vs. R.d: DEGs intersect in groups HC:PS vs. R and RIF:PS vs. R. Analysis of the frequency of differential ASEs demonstrated high occurrence rates across samples (Fig.  2 C). Evaluation of coding potential indicated that SE events most frequently retained coding capacity, whereas RI events led to frequent switches from coding to non-coding transcripts (Fig.  2 D). Gene Ontology enrichment analyses of ASEGs from both groups showed involvement in actin regulation, cell adhesion, GTPase activity, and RNA splicing (Fig.  2 E; Table S2 ). Corresponding DEGs were also enriched in cell adhesion and immune signaling pathways (Figure S2 B; Table S2 ). These processes have been implicated in embryo implantation 48 , especially adhesion junction proteins 49 , 50 and immune regulation 51 , 52 . Analysis of the PS and R phases in healthy controls identified 544 differential ASEs, including SE (n = 276, 50.73%), MXE (n = 127, 23.35%), A5SS (n = 35, 6.43%), A3SS (n = 52, 9.56%), and RI (n = 54, 9.93%). In the RIF: PS vs. R group, 804 ASEs were found, comprising 53.61% SE, 20.02% MXE, 6.59% A5SS, 9.70% A3SS, and 10.08% RI (Fig.  3 A). Fig. 3 Alternative splicing and immune cell analyses between PS and R phases in patients with RIF. (A) Analysis of differential ASEs between endometrial PS vs. R phases (blue) of HC and RIF patients (orange). The percentage of each type of ASEs within each comparison is represented by a pie chart (right) (SE, skipped exon; MXE, mutually exclusive exons; A3SS, alternative 3’ splice site; A5SS, alternative 5’ splice site; RI, retained intron). (B) Venn diagram showing the genes with alternative splicing events in HC: PS vs. R (cyan) and RIF: PS vs. R (blue). (C) The proportions of differential ASEs for each type in HC: PS vs. R (cyan) and RIF: PS vs. R (blue). (D) The frequency of possible changes in coding potential shown for SE, MXE, and RI, respectively. (E) Scatter plot displayed the GO analysis of genes with different alternative splicing events. BP: Biological process, CC: cellular component, MF: molecular function. (F) Comparing immune cell distribution between endometrial PS vs. R phase of RIF patients was shown in the boxplot (ns: p > 0.05, * p < = 0.05, ** p < = 0.01, *** p < = 0.001, **** p < = 0.0001). (G) Heatmap of the distribution of immune cells interacting with different ASEs in the endometrial R phases of RIF patients. Alternative splicing and immune cell analyses between PS and R phases in patients with RIF. (A) Analysis of differential ASEs between endometrial PS vs. R phases (blue) of HC and RIF patients (orange). The percentage of each type of ASEs within each comparison is represented by a pie chart (right) (SE, skipped exon; MXE, mutually exclusive exons; A3SS, alternative 3’ splice site; A5SS, alternative 5’ splice site; RI, retained intron). (B) Venn diagram showing the genes with alternative splicing events in HC: PS vs. R (cyan) and RIF: PS vs. R (blue). (C) The proportions of differential ASEs for each type in HC: PS vs. R (cyan) and RIF: PS vs. R (blue). (D) The frequency of possible changes in coding potential shown for SE, MXE, and RI, respectively. (E) Scatter plot displayed the GO analysis of genes with different alternative splicing events. BP: Biological process, CC: cellular component, MF: molecular function. (F) Comparing immune cell distribution between endometrial PS vs. R phase of RIF patients was shown in the boxplot (ns: p > 0.05, * p < = 0.05, ** p < = 0.01, *** p < = 0.001, **** p < = 0.0001). (G) Heatmap of the distribution of immune cells interacting with different ASEs in the endometrial R phases of RIF patients. Of the 342 ASEGs identified in the HC: PS vs. R group, only 5 (1.46%) overlapped with DEGs, whereas in the RIF: PS vs. R group, 38 (8.70%) of the 437 ASEGs were also DEGs (Figure S3 A). A total of 62 ASEGs were shared between these comparisons (Fig.  3 B), with 16 reported to be involved in embryo implantation (Table 3 ). Among these, APOL2, DUOXA1, SLC37A2, and ELP3 were also DEGs in both groups. AIMP1 and C3 have been associated with endometrial receptivity 7 , 46 . Frequency analysis showed high prevalence of ASEs across samples, with a trend towards higher frequency in RIF patients (Fig.  3 C). Coding potential analysis revealed that SE events generally preserved coding capacity, whereas RI events frequently resulted in non-coding isoforms (Fig.  3 D). GO enrichment analyses indicated that ASEGs were mainly related to cell adhesion, substrate junctions, and focal adhesion (Fig.  3 E; Table S3 ). DEGs were enriched in processes such as epidermis development, cytokine activity, and NK cell-mediated immunity (Figure S3 B; Table S3 ). DEGs identified in both healthy controls and RIF patients were enriched in immunomodulatory processes, including immune cell activation and migration. Analysis of immune cell composition revealed that CD8 + T cells, resting NK cells, and resting mast cells were significantly decreased, whereas monocytes and macrophages M0 were significantly increased in the receptive (R) phase endometrium of RIF patients (Fig.  2 F). Similarly, an increase in monocytes was observed in the R phase endometrium of healthy controls (Figure S2 C). Heatmap analysis demonstrated that several altered immune cell types, such as resting NK cells and resting mast cells, were closely associated with differential ASEs in both RIF and HC samples (Fig.  2 G; Figure S2 D). In the comparison between post-receptive (PS) and receptive (R) phases, immune cell analysis showed a significant increase in CD8 + T cells, M1 macrophages, and M2 macrophages, and a significant decrease in activated NK cells in RIF patients (Fig.  3 F). Although no statistically significant differences in immune cell proportions were found between PS and R phases in healthy controls, similar trends in CD8 + T cells and activated NK cells were observed (Figure S3 C). Among these immune cell types, activated NK cells displayed high correlation with specific ASEs during the receptive phase (Fig.  3 G; Figure S3 D). These findings suggest that dynamic changes in immune cell populations during the window of implantation may be closely linked to recurrent implantation failure. RNA-binding proteins (RBPs) are critical regulators of alternative splicing. Transcriptomic profiling identified seven RBP-related genes ( YBX2, CPEB2, IGF2BP2, IGF2BP3, KHDRBS3, BRUNOL5, and RBFOX1 ) that were differentially expressed in at least one comparison group (Fig.  4 A). Among these, CPEB2 has been implicated in trophoblast-related regulatory programs 53 , 54 . Notably, KHDRBS3 showed the highest degree of connectivity with differential ASEs in our correlation network (Fig.  4 C) and exhibited a consistent decrease in the receptive phase compared with the pre-receptive phase in both HC and RIF cohorts (Fig.  4 D), nominating it as a candidate splicing regulator associated with receptivity transitions. Fig. 4 The regulation network of ASEs by differentially expressed RBPs genes. (A) Differentially expressed RBP genes in four subgroups. (B) Significant enrichment of RBP motifs around AS events. On the x-axis, 5’ to 3’, the position (R) of each splicing event is shown relative to the event type. Each region is numbered from R1 to R(n), where n is the total number of regions for each event type. The binding motifs of the RBP are described on the y-axis. Each panel shows a single event type (SE, MXE, A5SS, A3SS, RI). The significance of the RBPs is indicated by the color and size of the circles. (C) Regulatory network of seven differentially expressed RBPs for the five types of alternative splicing (SE, MXE, A5SS, A3SS, RI) of ASEs. (D) Boxplots of KHDRBS3 expression in endometrial PR vs. R phase of HC population and RIF patients. p-values for each dataset are shown with corresponding group names. The regulation network of ASEs by differentially expressed RBPs genes. (A) Differentially expressed RBP genes in four subgroups. (B) Significant enrichment of RBP motifs around AS events. On the x-axis, 5’ to 3’, the position (R) of each splicing event is shown relative to the event type. Each region is numbered from R1 to R(n), where n is the total number of regions for each event type. The binding motifs of the RBP are described on the y-axis. Each panel shows a single event type (SE, MXE, A5SS, A3SS, RI). The significance of the RBPs is indicated by the color and size of the circles. (C) Regulatory network of seven differentially expressed RBPs for the five types of alternative splicing (SE, MXE, A5SS, A3SS, RI) of ASEs. (D) Boxplots of KHDRBS3 expression in endometrial PR vs. R phase of HC population and RIF patients. p-values for each dataset are shown with corresponding group names. Analysis of RBP binding motifs showed significant enrichment in regions flanking SE and MXE events (Fig.  4 B). Correlation analysis between expression levels of RBPs and ASE inclusion (PSI) values revealed that KHDRBS3 was associated with the highest number of ASEs, followed by IGF2BP2, YBX2, IGF2BP3, BRUNOL5, RBFOX1, and CPEB2 (Fig.  4 C). Notably, KHDRBS3 expression exhibited a consistent trend of significant reduction in the receptive phase compared to the pre-receptive phase in both HC and RIF groups (Fig.  4 D). These results support the hypothesis that RBPs contribute to splicing regulation in endometrial receptivity. Genes harboring pathogenic or aberrant splice variants are potential therapeutic targets. We collected 178 ASEGs shared between HC: PR vs. R and RIF: PR vs. R groups, and 62 ASEGs shared between HC: PS vs. R and RIF: PS vs. R groups, for drug repurposing analysis. Differentially spliced genes were uploaded to the Connectivity Map platform ( https://clue.io ), and the Drug Repurposing Hub 44 , part of the Connectivity Map, was used to identify small molecules with transcriptomic signatures inversely correlated with observed alterations. A signature-matching strategy was applied to retrieve compounds potentially capable of reversing disease-associated transcriptomic changes. Finally, the DrugBank database 45 was used to confirm compound approval status, pharmacological properties, and known indications. We identified 20 target ASEGs associated with 68 small molecules. Among them, GABRP, MAPK12 , and SGK1 were targeted by 26, 13, and 10 compounds, respectively (Table S4 ). Six ASEGs ( CD44, CKB, COX6C, NNMT, NR1H3, and ANXA1 ) previously reported in embryo implantation were matched with small molecules confirmed in DrugBank. These included hyaluronic acid ( CD44 ), creatine ( CKB ), cholic acid ( COX6C ), niacin ( NNMT ), T-0901317 ( NR1H3 ), and corticosteroids such as amcinonide, dexamethasone, and hydrocortisone phosphate ( ANXA1 ) (Fig.  5 A, B). Fig. 5 Drug repurposing analysis of genes with alternative splicing events. (A) Sankey diagram identifying small molecules targeting genes with differentially alternative splicing events. (B) Chemical structures of eight small molecule drugs. Drug repurposing analysis of genes with alternative splicing events. (A) Sankey diagram identifying small molecules targeting genes with differentially alternative splicing events. (B) Chemical structures of eight small molecule drugs. These compounds represent potential candidates for therapeutic intervention in RIF patients.

This prospective observational study was approved by the Ethics Committee of Shanghai Ji Ai Genetics and IVF Institute of Obstetrics and Gynecology, affiliated with Fudan University (JIAI E2019-04; JIAI E2020-015). All participants provided written informed consent prior to sample collection. All methods were carried out in accordance with relevant guidelines and regulations, including the Declaration of Helsinki and its later amendments. In total, 163 endometrial biopsy samples were included, comprising 90 samples from fertile healthy controls (HC) and 73 samples from patients with recurrent implantation failure (RIF). The HC cohort (n = 90) was obtained from our previously published endometrial receptivity transcriptome study, in which endometrial biopsies were collected across the natural-cycle luteinizing hormone (LH)–timed phases (LH + 3, LH + 5, LH + 7, and LH + 9) under standardized clinical procedures and subsequently profiled by RNA sequencing as reported previously 6 . These HC participants were fertile volunteers with proven fertility and no history of infertility, recurrent miscarriage, or known uterine pathology. The RIF cohort (n = 73) represents newly sequenced samples generated in the present study, collected from 35 patients undergoing a harmonized clinical management protocol at our center. Baseline demographic/clinical characteristics of the newly sequenced RIF cohort are summarized in Table 1 , and baseline hormone profiles under the HRT regimen are provided in Table S1 . Table 1 Basic characteristics of 73 samples from 35 RIF patients. Characteristic P + 3 (PR) P + 5 (R) P + 7 (PS) P value Endometrial samples (N) 24 26 23 Age (years): Mean ± SD 33.26 ± 4.51 33.48 ± 4.38 33.60 ± 4.39 NS BMI (kg/m 2 ): Mean ± SD 21.67 ± 2.17 21.76 ± 2.33 21.68 ± 2.46 NS Prior implantation failures (N): Mean ± SD 4.20 ± 1.81 4.11 ± 1.77 3.78 ± 1.20 NS Length of the menstrual cycle (Days): Mean ± SD 30.78 ± 2.66 32.00 ± 6.49 32.09 ± 6.78 NS Endometrium thickness (mm): Mean ± SD 8.95 ± 1.15 9.01 ± 1.47 9.15 ± 1.48 NS P + 3: 3rd day after starting progesterone administration; PR: pre-receptive; P + 5: 5th day after starting progesterone administration; R: receptive; P + 7: 7th day after starting progesterone administration;PS: post-receptive; BMI: Body Mass Index. The body mass index is the weight in kilograms divided by the square of the height in meters. Basic characteristics of 73 samples from 35 RIF patients. P + 3: 3rd day after starting progesterone administration; PR: pre-receptive; P + 5: 5th day after starting progesterone administration; R: receptive; P + 7: 7th day after starting progesterone administration;PS: post-receptive; BMI: Body Mass Index. The body mass index is the weight in kilograms divided by the square of the height in meters. RIF was defined as unexplained RIF with ≥ 3 failed embryo transfer attempts involving ≥ 4 high-quality embryos, where high-quality embryos were defined as either day-5 blastocysts graded at least 4BB (Gardner’s classification) or cleavage-stage embryos with ≥ 7 cells. Eligible patients were aged 20–40 years, had BMI 19–24 kg/m 2 , and an endometrial thickness of ≥ 7 mm at the time of endometrial preparation. To minimize confounding, patients with untreated hydrosalpinx, endometrial disease, severe adenomyosis, diminished ovarian reserve, genetic disorders, or immune abnormalities, or those unwilling to participate, were excluded from the study. For all RIF participants, embryo grading criteria and embryo transfer–related clinical procedures (including cycle management and luteal support strategy) followed uniform institutional standard operating procedures throughout the study period. For RIF patients, endometrial preparation was performed using artificial hormone replacement therapy (HRT). On menstrual cycle day 2, oral estradiol valerate (Progynova; Bayer, Leverkusen, Germany) was initiated at 4–6 mg/day until endometrial thickness reached ≥ 7 mm. Serum progesterone was then assessed; when progesterone was ≤ 1.5 ng/mL, 90 mg/day sustained-release vaginal progesterone gel (Crinone; Merck-Serono, Darmstadt, Germany) was administered. The start day of progesterone exposure was designated as P + 0, and endometrial biopsies were performed at P + 3, P + 5, and P + 7, corresponding to the pre-receptive, receptive, and post-receptive phases, respectively. Baseline gonadotropin and steroid hormone measurements associated with the HRT protocol are provided in Table S1 . Healthy controls underwent natural cycles with daily urinary LH monitoring (Eupregna, China). The day of the LH surge was defined as LH + 0, and biopsies were performed at LH + 3/LH + 5 (pre-receptive), LH + 7 (receptive), and LH + 9 (post-receptive), consistent with the timing scheme described in the prior publication 6 . Endometrial biopsy specimens were collected aseptically using a sterile, single-use endometrial suction catheter (Yudu Medical Apparatus and Instruments Co., Ltd., Suzhou, China). Immediately after collection, tissues were immersed in Allprotect Tissue Reagent (QIAGEN GmbH, Hilden, Germany), transported at 4 °C using standardized cold-chain procedures, and cryopreserved at –80 °C until downstream processing. Total RNA extraction, library preparation, and high-throughput sequencing were performed following the protocols described previously 6 . Sequencing reads were aligned to the GRCh38 human reference genome. Quality control, read mapping, and transcript quantification were conducted using established pipelines as previously reported 6 . All RNA-seq datasets generated in this study have been deposited in the Gene Expression Omnibus (GEO) under accession number GSE287072 . Quality-controlled data were analyzed using rMATS v4.1.1 to identify differential alternative splicing events (ASEs) across the following comparisons: (1) HC: PR vs. R; (2) HC: PS vs. R; (3) RIF: PR vs. R; and (4) RIF: PS vs. R. Five basic AS types were evaluated: skipped exons (SE), mutually exclusive exons (MXE), alternative 5′ splice sites (A5SS), alternative 3′ splice sites (A3SS), and retained introns (RI). The JCEC model (reads on target and junction counts) was applied for quantification. Significant ASEs were defined as those with an absolute inclusion level difference (ΔPSI) > 0.1 and false discovery rate (FDR) < 0.05, as recommended by Shen et al. 39 . Genes containing significant differential ASEs (ASEGs) were subjected to gene set enrichment analysis using the ClusterProfiler R package 40 . The MASER Bioconductor package was used to assess coding potential alterations resulting from splicing events 41 . Immune cell infiltration was estimated using CIBERSORT 42 , applying a validated signature matrix to deconvolute 22 immune cell types from transcriptomic data. Correlations between ASEs and immune cell proportions were assessed by Pearson correlation. The rMAPS2 platform was used to assess enrichment of RBP binding motifs around alternative splicing regions. Correlations between RBP expression and ASE inclusion levels were computed. Networks of RBPs and ASEs with correlation coefficient > 0.5 and p  < 0.05 were visualized using Cytoscape v3.9.1 43 . Differentially spliced genes were uploaded to the Connectivity Map platform ( https://clue.io ), and the Drug Repurposing Hub 44 , part of the Broad Institute’s Connectivity Map, was used to identify small molecules with transcriptomic signatures inversely correlated to the observed alterations. To prioritize already approved small molecule candidates, overlapping ASEGs from the HC: PR vs. R and RIF: PR vs. R groups were searched in the repurposing database. The same analyses were performed for overlapping ASEGs identified in the HC: PS vs. R and RIF: PS vs. R comparisons. A signature-matching strategy was applied to retrieve compounds potentially capable of reversing disease-associated transcriptomic changes. Finally, the DrugBank database ( https://go.drugbank.com ) 45 was used to confirm compound approval status, pharmacological properties, and known indications. Group differences were assessed using one-way analysis of variance (ANOVA). Results are presented as mean ± standard deviation (SD). A p -value < 0.05 was considered statistically significant. All analyses were performed using R (v4.1.3) and Python (v3.6.10).

This study demonstrates that endometrial receptivity is characterized by phase-dependent alterations in gene expression, alternative splicing events, and immune cell composition. The comprehensive transcriptomic profiling conducted here revealed that alternative splicing, particularly involving KHDRBS3 and other RNA-binding proteins, may play a pivotal role in regulating endometrial function during the window of implantation. Furthermore, the identification of differentially expressed genes and associated immune signatures highlights the complex molecular landscape underlying recurrent implantation failure. Collectively, these findings provide valuable insights into the mechanisms governing embryo implantation and may inform the development of novel diagnostic biomarkers and therapeutic strategies to improve implantation outcomes in assisted reproduction.

Recurrent implantation failure (RIF) remains a major challenge in assisted reproductive technology, and its pathogenesis is widely considered multifactorial, involving impaired endometrial receptivity, dysregulated immune–stromal crosstalk, and embryo‑related factors. In this study, we profiled bulk endometrial transcriptomes across receptivity phases to characterize coordinated changes in gene expression, alternative splicing (AS), and inferred immune composition in healthy controls (HC) and RIF patients. Across both cohorts, skipped exons (SE) and mutually exclusive exons (MXE) were the most prevalent AS types, consistent with the predominance of exon-skipping patterns reported in human tissues. We observed a larger shift in both differentially expressed genes (DEGs) and alternative splicing event genes (ASEGs) when comparing the pre-receptive (PR) to receptive (R) phase than when comparing the post-secretory (PS) to R phase, supporting extensive transcriptomic remodeling during establishment of the window of implantation. Functional enrichment of phase-associated ASEGs/DEGs implicated pathways linked to cell adhesion, actin cytoskeleton remodeling, and small-GTPase signaling—processes repeatedly connected to implantation competence through regulation of epithelial–stromal interactions, trophoblast attachment, and endometrial remodeling 48 , 50 , 51 . Together, these data support a model in which receptivity acquisition involves not only transcriptional reprogramming but also broad splicing rewiring affecting adhesion- and cytoskeleton-related networks. RNA-binding proteins (RBPs) are central regulators of splice-site choice, and their dysregulation can generate coordinated AS programs. By integrating RBP expression with AS patterns, we prioritized several RBPs as candidate upstream regulators of receptivity-associated splicing transitions. Among these, KHDRBS3 emerged as a notable candidate: it exhibited broad connectivity with differential AS events in the RBP–ASE network (Fig.  4 C) and showed a consistent decrease from PR to R in both HC and RIF cohorts (Fig.  4 D). KHDRBS3 has been implicated in RNA metabolism and signal-transduction-related processes in other biological contexts, suggesting that altered KHDRBS3 activity could plausibly reshape splicing programs relevant to endometrial remodeling. In our dataset, the timing of KHDRBS3 down-regulation coincided with the most pronounced ASEG/DEG transitions, supporting the hypothesis that KHDRBS3-centered splicing regulation participates in the molecular switch into the receptive state. Importantly, because bulk RNA-seq aggregates signals across multiple endometrial cell types, the observed KHDRBS3 dynamics may reflect both cell-intrinsic regulation and compositional shifts; cell type–resolved analyses will be required to localize the dominant cellular source(s) and to establish direct KHDRBS3 targets. Endometrial immune balance is another determinant of receptivity, and our CIBERSORT-based deconvolution suggested coordinated immune remodeling across phases. Notably, the inferred immune changes were phase-contrast dependent rather than uniform: in PR → R, resting NK cell signatures and several lymphocyte subsets decreased while monocytes and macrophage M0 increased in the receptive-phase endometrium of RIF patients (Fig.  2 F), whereas in PS → R, activated NK cell signatures decreased alongside shifts in macrophage subsets (Fig.  3 F). Across both comparisons, the most reproducible pattern involved the NK–myeloid axis, indicating an imbalanced remodeling trajectory during the transition toward receptivity in RIF. This pattern is biologically plausible because appropriate implantation requires synchronized immune tolerance and tissue remodeling, processes in which uterine NK cells and myeloid populations play key roles 30 – 32 , 53 . In parallel, we observed that altered immune fractions were correlated with specific differential ASEs in the receptive phase (Figs.  2 G and 3 G), supporting a link between splicing programs and immune microenvironment dynamics. To place the bulk immune findings into a more granular cellular context, recent single-cell endometrial atlases have mapped dynamic immune and stromal populations across the menstrual cycle, including resident NK and monocyte/macrophage lineages that vary with hormonal state 55 . In addition, single-cell studies focusing on RIF have reported perturbations in immune and stromal programs compared with fertile controls, including altered NK-cell states and macrophage/monocyte signatures in RIF endometrium 56 , 57 . Although single-cell studies differ in cohort characteristics and sampling windows, these reports provide a cell-type–resolved framework that is qualitatively consistent with the immune remodeling patterns inferred from our bulk deconvolution, while highlighting an important limitation of bulk RNA-seq/CIBERSORT: “NK cells” and “macrophages” at the bulk level likely reflect shifts in specific subsets and activation states rather than uniform lineage-wide changes. This is further supported by high-resolution single-cell mapping of the maternal–fetal interface, which has delineated functionally distinct NK subpopulations with divergent cytokine/chemokine and tissue-remodeling programs 58 . A salient implication of our integrated results is that AS regulation and immune remodeling may be mechanistically coupled. Splicing changes affecting adhesion molecules, cytoskeletal regulators, or immune signaling mediators could reshape cell–cell interactions and cytokine networks, thereby influencing recruitment or activation of NK cells and monocytes/macrophages. Conversely, inflammatory or stress cues derived from myeloid populations may affect RBP activity and spliceosomal regulation, reinforcing aberrant AS programs. In this framework, the KHDRBS3-associated splicing signature observed here may represent one molecular node at which receptivity-linked remodeling intersects with immune dysregulation in RIF. Beyond mechanistic inference, we applied a drug repurposing strategy to prioritize compounds targeting genes implicated in the AS/RBP network. Several candidate small molecules were identified, including agents with immunomodulatory or nuclear receptor-related activity, suggesting potential routes to restore receptivity-associated transcriptomic programs. Nevertheless, translating these computational candidates into clinical interventions will require careful preclinical validation, including assessment of endometrial bioavailability, timing relative to the window of implantation, and the directionality of immune effects. Several limitations warrant consideration. First, embryo aneuploidy cannot be fully excluded because preimplantation genetic testing was not performed; although RIF enrollment required repeated transfers of high-quality embryos, embryo-intrinsic factors remain a potential confounder. Second, bulk RNA-seq cannot definitively assign DEGs/ASEs (including KHDRBS3-linked splicing signals) to specific endometrial cell types, nor can it resolve NK and monocyte/macrophage subpopulations. Third, the RBP–AS associations are correlative; establishing causality and defining direct targets will require orthogonal validation (e.g., perturbation assays coupled with isoform-level readouts and/or RBP–RNA interaction assays). In summary, our data indicate that acquisition of endometrial receptivity is accompanied by marked phase-dependent remodeling of gene expression and alternative splicing, and that RIF is associated with coordinated perturbations in splicing regulatory RBPs (including KHDRBS3) and dysregulated immune remodeling involving the NK–monocyte/macrophage axis. These integrated observations refine the transcriptomic framework of implantation failure and support the development of receptivity-focused biomarkers and hypothesis-driven functional studies to delineate causal pathways linking splicing regulation, immune–stromal interactions, and implantation outcomes.

Although assisted reproductive technology (ART) has made significant progress over recent decades, 10% of patients still experience recurrent implantation failure (RIF), which remains a major challenge in clinical practice. RIF has been variably defined in the literature. According to the updated ESHRE guidelines (2023), a comprehensive assessment of contributing factors—including maternal age, embryo ploidy, and uterine pathology—is recommended when establishing a diagnosis of RIF. In the present study, we adhered to these recommendations and adopted an even more stringent definition: women under 40 years of age who failed to achieve a clinical pregnancy despite the transfer of at least four high-quality embryos across a minimum of three embryo transfer cycles 1 , 2 . As sequencing technologies evolve, transcriptomic data have become increasingly accessible. Studies have shown that the window of implantation (WOI), a short menstrual period when the endometrium is receptive to blastocyst transfer, is critical for successful implantation 3 , 4 . Transcriptomic analyses have revealed distinct genomic characteristics across the pre-receptive (PR), receptive (R), and post-receptive (PS) phases of the endometrium 5 , 6 . Based on these findings, transcriptome-based models to predict endometrial receptivity have been developed, leveraging transcript levels to offer personalized assessments of the WOI 6 – 10 . This approach has improved the diagnosis and management of RIF patients 11 , 12 . However, most existing studies have focused primarily on differential gene expression and have not comprehensively explored alternative splicing events (ASEs) as an additional regulatory layer of endometrial function. Alternative splicing (AS), which transforms a single mRNA precursor into multiple transcript variants, enhances proteomic diversity and regulates numerous cellular processes 13 , 14 . More than 60% of human genes contain multiple exons, many of which give rise to cell-type-specific isoforms 15 . Studies have demonstrated changes in AS during embryonic development and early pregnancy. By analyzing single-cell RNA-seq data, Tian et al. 16 found that AS is widespread in preimplantation embryo development, particularly at the two-cell stage. A homozygous splicing mutation of HFM1 demonstrated the potential risk of chromosomal anomalies under the RIF phenotype 17 . Similarly, in human embryonic arrest and RIF patients, splicing mutations have been shown to cause abnormal alternative splicing resulting in truncated MEI1 proteins 18 . Furthermore, Zhao et al. found that increased HOXA11-AS expression caused impaired PKM2 splicing and attenuated decidualization, a change consistently observed in RIF patients 19 . RNA-binding proteins (RBPs) are a diverse protein family capable of binding single- or double-stranded RNA 20 – 22 . RBPs can promote the inclusion or skipping of exons by binding to splicing regulatory elements. Changes in the levels and activity of RBPs can cause dysregulation of AS 23 , 24 . Several studies have reported that RBPs play an important role in embryo implantation through the regulation of AS events, particularly in endometrial decidualization 19 , 25 – 28 . For example, PTBP1 has been shown to regulate PKM1/2 alternative splicing and gene expression, thereby affecting decidualization in RIF patients 19 . Further analysis of transcriptome data may be essential to understand the distribution, complexity, and regulation of alternative splicing across endometrial receptivity phases and to identify candidate RBPs involved in the control of key splicing events in endometrial function. Endometrial immune dysfunction is another factor reducing receptivity and contributing to implantation failure. Various immune cells in the endometrium, such as natural killer (NK) cells, macrophages, and T cells, are essential for regulating receptivity and embryo implantation 29 , 30 . Uterine NK cells, T regulatory (Treg) cells, dendritic cells (DCs), and macrophages can directly or indirectly influence uterine epithelial adhesion, stromal cell transformation, trophoblast differentiation and invasion, and uterine vascular adaptation 31 – 33 . These immune cell profiles have been used to assess endometrial receptivity 34 – 36 . Although earlier studies have suggested that intravenous immunoglobulin therapy might benefit women with RIF 37 , 38 , recent ESHRE guidelines have questioned this practice, and its efficacy remains controversial. Therefore, the connection between immune cell composition and endometrial receptivity warrants further investigation with updated clinical evidence. To investigate transcriptome dynamics during embryo implantation, we collected endometrial tissue at pre-receptive, receptive, and post-receptive phases from healthy controls (HC) with proven fertility and patients with RIF. RNA-seq was performed on all samples, followed by an extensive transcriptomic analysis, including DEGs, ASEs, immune cell estimation, and RBP network construction, with the aim of elucidating their regulatory relationships in endometrial function. By comprehensively characterizing the endometrial transcriptome across receptivity phases, this study aimed to clarify the biological significance of widespread alternative splicing and further evaluate immune cell changes and their potential role in implantation failure. In addition, an ASE–RBP correlation network was constructed to explore the dysregulation of splicing by RBPs, and drug repurposing analysis was performed to identify potential therapeutic candidates targeting dysregulated AS. Compared to previous transcriptomic studies, we conducted further analysis of ASEs in endometrial receptivity and attempted to identify key RBPs responsible for splicing regulation. By revealing transcriptional diversity and dysregulation, we aimed to provide a valuable resource for studying splicing and identifying potential biomarkers and therapeutic targets in RIF, as summarized in the study design flowchart (Fig.  1 ). Fig. 1 Overview of the study design and analytical workflow. Endometrial tissue was collected from pre-receptive (PR), receptive (R), and post-receptive (PS) phases from healthy fertile controls (HC) and patients with recurrent implantation failure (RIF). Bioinformatic analysis included differential gene expression (DEG) and alternative splicing event (ASE) identification, immune cell estimation, RNA-binding protein (RBP) motif enrichment, and drug repurposing analysis. The figure was created with BioRender.com. Overview of the study design and analytical workflow. Endometrial tissue was collected from pre-receptive (PR), receptive (R), and post-receptive (PS) phases from healthy fertile controls (HC) and patients with recurrent implantation failure (RIF). Bioinformatic analysis included differential gene expression (DEG) and alternative splicing event (ASE) identification, immune cell estimation, RNA-binding protein (RBP) motif enrichment, and drug repurposing analysis. The figure was created with BioRender.com.

Supplementary Information 1. Supplementary Information 2. Supplementary Information 3. Supplementary Information 4. Supplementary Information 5. Supplementary Information 6. Supplementary Information 7. Supplementary Information 1. Supplementary Information 2. Supplementary Information 3. Supplementary Information 4. Supplementary Information 5. Supplementary Information 6. Supplementary Information 7.