Comparing transcriptome profiles of human embryo cultured in closed and standard incubators

Ethics statement

This study was approved by the Institutional Review Board (IRB) of Chongqing Health Center for Women and Children (2016-RGI-01). We followed the guiding principles from the Ministry of Science and Technology (MOST) for the review and approval of human genetic resources. All donor couples voluntarily donated embryos after signing written informed consent at the Chongqing Reproductive and Genetics Institute in the Chongqing Health Center for Women and Children.

Intracytoplasmic sperm injection

We performed Intracytoplasmic sperm injection (ICSI) within 5 h of oocyte retrieval under an inverted microscope (Olympus IX70; Olympus Optical Co. Ltd., Tokyo, Japan) with a micromanipulation system (CellTram® 4r; Eppendorff, Hamburg, Germany). During ICSI, we placed the oocytes in pre-equilibrated culture droplets and covered them with six ml of mineral oil (Ovoil, Vitrolife, Gothenburg, Sweden). We used a holding pipette to position the oocyte. The sperm was injected into the cytoplasm with a micropipette when the first polar body reached the 6 or 12 o’clock position.

Embryo culture

For the SI group, we placed the injected oocytes into a pre-equilibrated culture dish (Thermo Scientific, Waltham, MA, USA) with 25 µL of culture droplets (Vitrolife Sweden AB, Gothenburg, Sweden) covered with 1.2 mL of paraffin oil (Vitrolife Sweden AB, Gothenburg, Sweden). The embryos were cultured in an SI (MCO-5M; Sanyo, Osaka, Japan) at 37 °C with 5% O₂ and 6% CO₂ until embryo transfer on Day 3.

For the TLI group, we transferred the injected oocytes into pre-equilibrated Embryoslides and cultured them in the EmbryoScope™ until the time of transfer. The embryos cultured in the EmbryoScope™ had the same culture conditions as the SI group.

Embryo scoring

We performed morphological assessments of the SI- and TLI-cultured embryos at the same time points using the same criteria. Additional information obtained during TLI was not used for embryo assessment or selection. In this study, the time points for morphological assessments followed the guidelines suggested by the ESHRE/Alpha consensus (31). The embryos were assessed at 17 ± 1, 44 ± 1, and 68 ± 1 h after ICSI. We scored embryos on their blastomere shape, blastomere number, and fragmentation rate. An embryo was defined as grade 1 when it had an even blastomere shape and <10% fragmentation, grade two when it had uneven blastomeres and 10–25% fragmentation, grade three when it had uneven blastomeres and 25–35% fragmentation, and not recommended for transfer when the fragmentation >35%. Generally, if one patient have more than eight embryos with grade 1–3, we will transfer two embryos with top quality, and freeze the other top six embryos. Then, the other embryos excluding the eight embryos were defined as not reach the transfer degree in this cycle, even if embryos with grade three. However, if one patient did not have enough embryos with grade one or two, for example, only one embryo with grade three. We will only transfer one embryo with grade three. Here, the embryos used in this study were all from the first situation, and were grade three with eight blastomere.

Collecting human zygotes and eight-cell embryos

The patients who donated embryos were between 25 and 30 years old without a history of genetic diseases or smoking. We collected zygotes 6 h after ICSI, and eight-cell embryos on Day 3. Each couple donated only one zygote or embryo. To avoid possible embryo heterogeneity, each RNA-Seq pool had five zygotes or embryos from five different patients. We processed a total of 15 zygotes and 30 embryos (15 TLI-cultured and 15 SI-cultured) for RNA-Seq. All embryos were frozen using vitrification and were then stored. The embryos were cultured for 2 h after thawing, briefly exposed to an acidic PBS solution for 5–10 s to remove the zona pellucida, and were washed three times in PBS. After all the embryos in one group were prepared, we removed them from the PBS, and immediately placed five embryos per tube into lysis buffer. All embryos used in this study exhibited normal morphology without reaching the transfer degree in this cycle.

RNA sequencing library generation

We performed amplification using the Smart-seq2 method. We used the Qubit® 3.0 Fluorometer and Agilent 2100 Bioanalyzer to check the quality of the cDNA product and to ensure its length was around 1–2 kb. The library was prepared following the manufacturer’s instructions (Illumina. Cat. FC-131-1024). After library preparation, we checked the library quality using the PerkinElmer LabChip® GX Touch and Step OnePlus™ Real-Time PCR System. The libraries were then sequenced on the Illumina HiSeq 4000 with a 150-bp paired-end.

RNA sequence data processing

We used Trim_Galore to remove raw sequence reads that contained adapter contamination and poor-quality reads with low PHRED scores. Tophat (version 2.4.0) was used to map all data to the hg38 genome. We used RSEM (version 1.3.0) to calculate the gene expression levels based on the GRCh38.92 annotation database. We performed differential expression analysis using DESeq2, and only kept genes with significantly higher changes (log2 fold changes ≥1 and false discovery rates (FDR) <0.05). We deposited all RNA-Seq data in the NCBI Gene Expression Omnibus (GEO) under accession number GSE114771.

Gene ontology analysis

We performed gene ontology (GO) analysis using the DAVID web-tool. Biological process, molecular function, and cellular component terms were selected as the background gene sets.

RT-qPCR

Total RNA extraction, reverse transcription, and pre-amplification were performed using the Single Cell-to-CT qRT-PCR Kit (4458237; Invitrogen, Carlsbad, CA, USA) according to the manufacturer’s instructions. We conducted real time PCR with SYBR Premix Ex Taq™ (TaKaRa, Kyoto, Japan) and the CFX96 Touch Real-Time PCR System (Bio-Rad, Hercules, CA, USA). We measured Ct values using Sequence Detection System software (Bio-Rad, Hercules, CA, USA), and analyzed the gene expression level using the 2^−ΔΔCt method. All experiments contained at least three biological replicates. We used of the average of GAPDH and HPRT for normalizing RT-qPCR data. The primers used in this study are shown in Table S1.

Statistical analysis

We used SPSS statistical software 19.0 for the statistical analyses, and one-way ANOVA to assess the differences between the two groups. For all analyses, differences were considered significant when p < 0.05.

Deep sequencing of the human embryo transcriptome

Using RNA-Seq technology, we studied the transcriptomic profiles of the three groups of human embryos: zygotes, SI-cultured, and TLI-cultured eight-cell embryos. After removing reads with adaptors, reads containing poly N, and low-quality reads from the raw data, a total of 68,950,172 (88.1% of the raw data), 32,041,678 (72.35%), and 27,625,662 (62.63%) trimmed reads remained in the zygote group; 65,996,908 (88.94%), 41,141,300 (70.82%), and 52,110,290 (83.02%) remained in the SI group; and 53,242,892 (80.21%), 45,191,346 (81.37%), and 39,410,428 (81.86%) remained in the TLI group, respectively (Table 1).

We analyzed the trimmed reads via TopHat for alignment with the reference genome (hg38). We mapped approximately 400 million paired-end reads derived from the nine samples, with a mapping rate ranging from 89% to 93% (Table 1). The typical number of detectable genes in the individual samples ranged from approximately 18,000 to 27,000 genes.

Global transcriptomic profile analysis

The number of trimmed reads mapped to a gene was normalized to transcripts per kilobase million (TPM) values to accurately measure the gene expression level. To determine whether the TLI- and SI-cultured embryos had distinct transcriptomic landscapes, we analyzed the RNA-Seq data using unsupervised hierarchical clustering. We included all genes that were expressed in at least one of the nine samples with TPM values ≥0.1. The results showed that embryos in the same developmental stage clustered together. Moreover, the TLI samples were similar to the SI samples and were distinct from the zygotes (Fig. 1A). Our principal component analysis (PCA) also found these distinct expression patterns (Fig. 1B).

Differentially expressed gene analysis

Although the global transcriptomic profiles of the TLI group were more similar to those of the SI group than the zygotes, we could not distinguish any small molecular differences between TLI- and SI-cultured embryos. Therefore, we independently compared the TLI and SI groups and identified a total of 539 differentially expressed genes (DEGs, FDR <0.05, >2-fold change). We verified the accuracy using cluster analysis, which successfully distinguished the TLI from the SI group (Fig. 2A). Of the 539 DEGs, 24 genes were exclusively expressed in the TLI group, and 39 genes were exclusively expressed in the SI group. Of the remaining 476 genes, 188 showed higher expression in the TLI group while 288 genes showed higher expression in the SI group (Table S2). We also detected the differential expression of solute carrier family (SLC) genes, including SLC5A12, SLC9B2, and SLC35A3, which play important roles in early embryonic development (Song et al., 2005). Figure 2B shows the fold-change breakdown across the DEGs of the two groups.

Functional analysis

To identify the functions of DEGs, we performed GO enrichment analysis of SI- and TLI-overexpressed genes. In SI-overexpressed genes, we found many more enriched GO terms, particularly biological process terms (83%, Fig. 3). The overexpressed genes in SI tended to execute functions, such as transcription, RNA splicing, and DNA repair. However, the overexpressed genes in TLI belonged to several different categories, including the cell differentiation and methyltransferase activity pathways. In spite of the similarities between the SI and TLI groups, there were still small differences that could influence specific pathways and molecular functions during embryonic development.

Validating differentially expressed genes

To verify the reliability of our sequencing data, we selected four upregulated and downregulated genes from several enriched pathways, including RNA splicing, DNA repair, cell differentiation, and methyltransferase activity. We used qPCR to analyze their expression levels in eight-cell embryos cultured in SI and TLI. We detected the upregulated expressions of RBM14, TPRAP, SMARCB1, and ELAVL1 in SI, and ANGPT1, SLC35A3, PRDM5, and ZIC3 in TLI (Fig. 4). These results were consistent with our sequencing data.