Implications
Over the last 40 to 50 years there has been a significant decline in reproductive performance of dairy cattle, due in part to low fertilization rates and early embryonic loss. Although genetics plays an important role in embryo development, the identification of genes associated with embryo quality and pregnancy success has been a serious challenge. In this study, we characterized several genes that showed altered expression in preimplantation embryos that do not develop normally to the blastocyst stage. These genes can be used as candidate markers for embryo quality in cattle.
Introduction
In the last few decades, production of embryos by in vitro fertilization (IVF) has become a common practice for the treatment of infertility in both humans and cattle, and although IVF procedures are considered to be a true revolution in the field of fertility, the success rate of IVF pregnancies is still low. For example, in humans, ∼8 out of 10 transferred embryos fail to implant successfully (Bromer and Seli, Reference Bromer and Seli2008), most likely because of the inability to identify embryos with the best development potential. In cattle, there has been a significant decline in fertility over the past few decades (Dobson et al., Reference Dobson, Smith, Royal, Knight and Sheldon2007; Leroy et al., Reference Leroy, Opsomer, Van Soom, Goovaerts and Bols2008) because of low fertilization rates and embryonic loss among other factors (Santos et al., Reference Santos, Thatcher, Chebel, Cerri and Galvão2004; Morris and Diskin, Reference Morris and Diskin2008). In order to identify genes affecting embryonic loss, and hence infertility, an IVF experimental system in cattle has been established in our laboratory and several candidate genes have been reported to be associated with embryo development (Khatib et al., Reference Khatib, Huang, Wang, Tran, Bindrim, Schutzkus, Monson and Yandell2009; Huang and Khatib, Reference Huang and Khatib2010; Huang et al., Reference Huang, Kirkpatrick, Rosa and Khatib2010a and Reference Huang, Yandell and Khatib2010b).
The recent high-throughput sequencing technologies, such as next-generation RNA sequencing (RNA-Seq), allow for characterization of the expression activities of thousands of transcripts at unprecedented resolution. In addition to nucleotide sequencing, with RNA-Seq one has the ability to estimate expression levels of transcripts, identify unannotated genomic regions and detect alternative splicing events. Recently, we reported the first application of RNA-Seq technology in cattle for the characterization of transcriptomic landscapes of embryos at the blastocyst stage (Huang and Khatib, Reference Huang and Khatib2010). Sequencing data in fastq format can be accessed by Gene Expression Omnibus (GEO) with the accession number GSE25082. We have identified 1 648 307 transcribed bases that have not been annotated in the bovine genome. Interestingly, of these novel transcriptional units, 25 transcripts were found to be differentially expressed between cattle embryos with different morphology and developmental statuses (degenerate embryos and normally developed blastocysts). However, it is well recognized that high-throughput analysis studies are hypothesis generating and that results of these studies should be validated in subsequent experiments. Thus, the first objective of the current study was to characterize these novel transcripts and validate their association with early embryo development using quantitative real-time PCR (qRT-PCR).
Recent studies have established that alternative splicing is widespread and plays a significant role in shaping phenotypic complexity (Blencowe, Reference Blencowe2006). Wang et al. (Reference Wang, Sandberg, Luo, Khrebtukova, Zhang, Mayr, Kingsmore, Schroth and Burge2008) estimated that more than 92% of human genes undergo alternative splicing and that most events of alternative splicing are regulated differently among tissues. The extent of alternative splicing in cattle genes is not known, and only about 83 000 exons have been characterized compared with about 200 000 exons in human genes. In a previous study, we utilized RNA-Seq data to characterize patterns of alternative splicing in blastocysts and degenerate embryos. We detected thousands of alternative splicing events including ‘exon skipping’, ‘alternative 5′ splice site’ and ‘alternative 3′ splice site’ in multi-exon genes and concluded that alternative splicing is widespread in bovine embryos. We found 23 splicing events that were significantly associated with embryonic developmental status.
For better understanding of the contribution of alternative splicing to early embryonic development, there is a need to validate and further investigate differential expression results obtained in the high-throughput analysis (Huang and Khatib, Reference Huang and Khatib2010) at the single-gene level and to biologically replicate these results. Thus, the second objective of this study was to experimentally estimate the differential expression of the splice variants of a subset of genes undergoing alternative splicing in cattle embryos. We hypothesize that aberrant expression of genes in the developing embryo could lead to abnormalities such as growth retardation and mortality, thereby affecting the reproductive performance of the parent individual.
Material and methods
Characterization of novel bovine sequences and alternatively spliced genes
In a previous study, we utilized the RNA-Seq technology to identify novel transcripts – those that have not previously been annotated – in two groups of embryos (normal blastocysts and degenerate embryos) that differ in their development and morphology (Huang and Khatib, Reference Huang and Khatib2010). The criteria, definitions and discovery of the novel transcribed units (TUs) were described in Huang and Khatib (Reference Huang and Khatib2010). In brief, they were defined as regions containing transcripts that were at least 1 kb away from known gene boundaries, of length more than 500 bp, and covered by at least 40 fragments per kb. A total of 1 648 307 transcribed unannotated bases were identified in which 1272 and 1348 novel TUs had met these criteria in blastocysts and degenerate embryos, respectively. A total of 25 novel TUs were found to be differentially expressed between blastocysts and degenerate embryos. Of these, eight TUs were chosen for this study for further investigation and validation of RNA-Seq results. The criteria for selecting these sequences included sequence similarity with known human genes, high counts in RNA-Seq analysis and close vicinity to known cattle genes. The basic local alignment search tool (BLAST) function was used to search for sequence similarities between sequences of the novel transcripts and sequences of human and bovine (BTAU 4.0) deposited in National Center for Biotechnology Information (NCBI) web site (www.ncbi.nlm.nih.gov).
For validation and further characterization of alternative splicing results obtained in RNA-Seq analysis (Huang and Khatib, Reference Huang and Khatib2010), the genes MYL6, NOP10, RNF187, RPS24 and RPS28 were chosen because they each had more than one splice variant and because they showed differential alternative splicing between blastocysts and degenerate embryos. Figure 1 shows the genomic structure of the alternative transcripts and the locations of primers used in the expression analysis.
Genomic structure of splice variants used in the expression analysis. F and R denote the positions of the forward and reverse primers, respectively, used in the quantitative real-time PCR reactions. START and STOP mark the positions at which translation into protein sequences begins and ends, respectively.
IVF
For the validation of expression of novel TUs and alternatively spliced genes obtained in the RNA-Seq analysis, two types of embryos that differed in morphological and developmental phenotypes were produced as described in Khatib et al. (Reference Khatib, Monson, Schutzkus, Kohl, Rosa and Rutledge2008) and in Huang and Khatib (Reference Huang and Khatib2010). Embryos that showed a distinct inner cell mass and blastocoele were classified as ‘blastocysts’, whereas those that failed to make transition from morula to blastocyst stage were classified as ‘degenerates’. Three biological replicates of pairs of embryo pools were collected; each consisted of 20 blastocysts or degenerate embryos.
qRT-PCR
Total RNA was extracted from pools of embryos using the RNaqueous Micro kit (Ambion, Austin, TX, USA) and quality controlled using a RNA6000 PicoChip (Agilent Technologies, Santa Clara, CA, USA). Primers (Table 1) were designed to amplify fragments in the qRT-PCR reactions using the Beacon Designer software (Premier Biosoft International, Palo Alto, CA, USA) and by NCBI Primer-BLAST for specific amplification. Primers were designed to span exon–exon boundaries to avoid co-amplification of genomic DNA.
Primer sequences and amplicon sizes of the novel sequences and alternatively spliced genes used in the qRT-PCR
qRT-PCR = quantitative real-time PCR; GAPDH = glyceraldehyde-3-phosphate dehydrogenase.
An internal control gene for normalization was chosen as described in Vandesompele et al. (Reference Vandesompele, De Preter, Pattyn, Poppe, Van Roy, De Paepe and Speleman2002). Briefly, this method assumes that the ratio of two internal control genes is uniform across samples, regardless of the cell type or environmental condition. The control gene producing the smallest relative stability value M (the average of the pair-wise variation when compared with the other control genes) would show the most stable expression across samples, assuming no co-regulation and is thus chosen for normalization. The housekeeping genes beta-actin (ACTB), acidic ribosomal protein large (RPLP0) and glyceraldehyde-3-phosphate dehydrogenase (GAPDH) were tested across blastocyst and degenerate samples with GAPDH having the lowest calculated stability value (M = 0.26). Reactions of qRT-PCR were carried out in a DNA Engine-Opticon 2 Detection System (MJ Research, Watertown, MA, USA) using the iQ SYBR Green Supermix kit (Bio-Rad Laboratories, Hercules, CA, USA). All samples were run in quadruplicates. Relative gene expression values were calculated using the 2−ΔΔCt method (Livak and Schmittigen, Reference Livak and Schmittgen2001). For expression analysis, normalized gene expression values (ΔCt) were analyzed using a general linear model including the fixed effects of the bull, the type of embryo (blastocyst or degenerate embryos) and the random effect of the pool. The mean and the range of the fold change for each gene were calculated as 2−ΔΔCt using the estimated ΔΔCt value ± s.e.
Results and discussion
Differential expression of novel transcripts in bovine embryos
Table 2 shows the differential expression of the novel transcripts in blastocysts and degenerate embryos and their sequence homology to known human and cattle genes. Although the type of embryo (blastocyst or degenerate embryos) showed significant association with gene expression in each of the analyses, we did not detect significant association between the bull effect and differences in gene expression in any of the analyses. A range of 2.5- to 90-fold difference in expression was observed between embryo groups; all the novel transcripts were found to be upregulated in degenerate embryos compared with blastocysts. Supplementary Table 1 shows results of BLAST analysis of novel transcripts including transcript length, coverage and sequence identity.
Characteristics of novel transcripts and their fold difference in degenerate embryos v. blastocysts
The transcript XLOC_000727 showed a 90-fold higher expression in degenerate embryos v. blastocysts (Table 2). A search in the NCBI database for sequences similar to XLOC_000727 revealed high similarity (78%) to the human zinc finger imprinted 3 (ZIM3) gene. ZIM3 is located in the PEG3 imprinted domain and comparative analysis with species including human, mice and cow indicated that this domain structure has been well preserved throughout mammalian evolution (Kim et al., Reference Kim, Bergmann, Choo and Stubbs2007). The mouse Zim3 is imprinted with maternal expression in adult testes and based on the lack of open reading frame, it was assumed that the gene has no protein-coding capability (Kim et al., Reference Kim, Bergmann, Wehri, Lu and Stubbs2001). In contrast, the human and bovine ZIM3 is protein-coding gene, but the imprinted status has not been determined (Kim et al., Reference Kim, Bergmann, Choo and Stubbs2007).
Despite the fact that the biological role of ZIM3 is unknown, it is well documented that duplications or mutations in the PEG3 domain result in prenatal lethality in mammals (Kuroiwa et al., Reference Kuroiwa, Kaneko-Ishino, Kagitani, Kohda, Li, Tada, Suzuki, Yokoyama, Shiroishi, Wakana, Barton, Ishino and Surani1996; Flisikowski et al., Reference Flisikowski, Venhoranta, Nowacka-Woszuk, McKay, Flyckt, Taponen, Schnabel, Schwarzenbacher, Szczerbal, Lohi, Fries, Taylor, Switonski and Andersson2010). This suggests a regulatory mechanism affecting prenatal development through PEG3 expression and therefore a possible relationship between ZIM3 and embryonic development. The upregulation of ZIM3 in degenerate embryos observed in this study is of particular interest given that imprinted genes have important roles in embryonic, placental and neonatal growth (Bressan et al., Reference Bressan, De Bem, Perecin, Lopes, Ambrosio, Meirelles and Miglino2009) and that aberrant expression of these genes has been reported to affect early embryonic development. For example, Tunster et al. (Reference Tunster, Tycko and John2010) showed that a doubling of expression of the imprinted gene Phlda2 in mice resulted in a significant reduction in embryonic glycogen stores and embryonic weight.
Using the BLAST function, XLOC_001837 transcript showed an 87% sequence similarity to isoform 8 of the sorbin and SH3 domain containing 2 (SORBS2) human gene, but no significant similarity to the bovine gene was found. Indeed, an Ensembl search revealed 65 splice variants of the human gene, whereas only two variants were found for the bovine gene (http://uswest.ensembl.org/), an indication of the incompleteness of the annotation of the bovine genome. The transcript XLOC_001837 was found to be upregulated in degenerate embryos (12-fold) compared with blastocysts (Table 2). This result is consistent with a previous report in which SORBS2 was found to be upregulated by the preimplantation factor in early stages of embryonic development (Paidas et al., Reference Paidas, Krikun, Huang, Jones, Romano, Annunziato and Barnea2010).
The transcript XLOC_001206 located 3 kb downstream of the SOX9-like bovine gene showed a 38-fold difference in expression with more of this transcript present in degenerate embryos compared with blastocysts (Table 2). This novel transcript showed 86% sequence similarity to the human sex determining region Y (SRY)-box 9 (SOX9) gene. SOX9 is a member of the SOX transcription factors family that shares the high mobility group DNA-binding motif with the mammalian testis-determining factor SRY. The gene has been reported to have a wide range of functions in differentiation and development including chondrogenesis and skeletal development, neural crest differentiation, hair follicle function, pancreas development and sex determination (Lefebvre et al., Reference Lefebvre, Behringer and de Crombrugghe2001). Furthermore, recent studies have provided accumulating evidence of a role for SOX9 in early embryonic growth. For example, downregulation of SOX9 expression by microinjection of siRNA resulted in small-sized embryos, abnormal development and early death of mouse embryos, which suggests that the gene is a regulator of embryonic growth during early development (Grandjean et al., Reference Grandjean, Gounon, Wagner, Martin, Wagner, Bernex, Cuzin and Rassoulzadegan2009). Thus, the differential expression of XLOC_001206 transcript observed in this study may provide further support for the role of SOX9 in early embryonic development.
The transcript XLOC_001203 was upregulated in degenerate embryos v. blastocysts (Table 2). BLAST function revealed 100% sequence similarity between this transcript and intron 7 of the bovine matrix-remodelling associated 7 (MXRA7) gene. The specific functions of MXRA7 in cattle are not known. However, in humans the gene was reported to be upregulated in the ectopic endometrium of women with endometriosis (Meola et al., Reference Meola, Rosa e Silva, Dentillo, da Silva, Veiga-Castelli, Bernardes, Ferriani, de Paz, Giuliatti and Martelli2010), and a gene ontology search revealed possible functions in tissue remodeling (Veiga-Castelli et al., Reference Veiga-Castelli, Silva, Meola, Ferriani, Yoshimoto, Santos, Squire and Martelli2010).
The transcript XLOC_000667, which showed an 86-fold higher level of expression in degenerate embryos v. blastocysts, had some sequence similarity with the F-box protein-like bovine gene (Table 2). The F-box gene family encodes proteins that contain at least one F-box domain, a protein structural motif of about 50 amino acids that mediates protein–protein interactions. The F-box proteins are involved in the ubiquitin system of protein degradation (Craig and Tyers, Reference Craig and Tyers1999). Moreover, these proteins have been associated with cellular functions such as signal transduction and regulation of the cell cycle (Kipreos and Pagano, Reference Kipreos and Pagano2000). Thus, the upregulation of transcript XLOC_000667 in degenerate embryos could indicate some function for this transcript in the cell cycle and the degradation of proteins of the developing embryo.
The transcripts XLOC_001523, XLOC_001513 and XLOC_000024 were upregulated in degenerate embryos compared with blastocysts (Table 2). A BLAST search revealed that these sequences have no similarity to known human, mouse or bovine sequences. Given that the mouse and the human genomes are relatively well annotated, these transcripts are most likely bovine specific. Transcript XLOC_000024 is located 74 kb upstream of the zinc finger protein 131 (ZNF131) gene, but has no homology to known genes. Recently, there is accumulating evidence that ∼93% of the human genome is transcribed into primary transcripts, which include a large number of unannotated transcripts without coding capacity (The ENCODE Project Consortium, 2007). Furthermore, it has been proposed that this pervasive transcription appears to be a universal characteristic for higher eukaryotic genomes (Dinger et al., Reference Dinger, Amaral, Mercer and Mattick2009). Thus, some of the novel sequences identified in this study could be transcribed from nongenic regions or from noncoding 5′ or 3′ untranslated regions of coding genes, which in turn play significant biological roles.
Differential expression of alternatively spliced transcripts
Figure 2 shows the expression profiles of splice variants of genes found to be differentially expressed between blastocysts and degenerate embryos. Transcripts 660 and 691 of the myosin light chain 6 (MYL6) gene showed 2.5- and 8.1-fold difference in expression, respectively, in degenerate embryos compared with blastocysts (Figure 2). These results are consistent with previous reports on the association of MYL6 with blastocyst formation. Goossens et al. (Reference Goossens, Van Soom, Van Poucke, Vandaele, Vandesompele, Van Zeveren and Peelman2007) reported that MYL6 is highly expressed at the blastocyst stage compared with earlier stages of development and concluded that this gene could be a candidate marker for blastocyst formation. However, no data has been reported on the involvement of the different splice variants of MYL6 in early embryonic development. The two alternative splicing isoforms of MYL6 differ by eight residues at the C terminal of the polypeptide chain and this pattern of alternative splicing is conserved between the human and the bovine (Huang and Khatib, Reference Huang and Khatib2010).
Fold difference in expression of alternative splicing transcripts in degenerate embryos v. blastocysts. Error bars indicate the range (minimum–maximum) of the fold difference.
Both transcripts 2939 and 3011 of the NOP10 ribonucleoprotein homolog (yeast) gene (NOP10) showed 11.2- and 5.7-fold higher expression, respectively, in degenerate embryos compared with blastocysts. NOP10 is a member of the H/ACA snoRNP gene family that is involved in the processes of ribosome biogenesis and telomere maintenance. A mutation in human NOP10 is associated with an autosomal recessive form of dyskeratosis congenita (Walne et al., Reference Walne, Vulliamy, Marrone, Beswick, Kirwan, Masunari, Al-Qurashi, Aljurf and Dokal2007). Furthermore, telomerase activity is dynamically regulated in early bovine embryonic development and appears to reach a peak at the blastocyst stage (Xu and Yang, Reference Xu and Yang2000). Thus, results of this study provide further evidence of the association of NOP10 with early embryonic development.
The splice variants of the genes ring finger protein 187 (RNF187), ribosomal protein s24 (RPS24) and ribosomal protein s28 (RPS28) also showed significant differences in expression between blastocysts and degenerate embryos. Interestingly, the functions of these genes are largely unknown in cattle. RPS24 and RPS28 are both components of the small 40S subunit of ribosome. Apart from their structural role in forming ribosomes, there are no reports on other functions of these genes. Nevertheless, mutations in human RPS24 were reported to be associated with Diamond–Blackfan anemia (Gazda et al., Reference Gazda, Grabowska, Merida-Long, Latawiec, Schneider, Lipton, Vlachos, Atsidaftos, Ball, Orfali, Niewiadomska, Da Costa, Tchernia, Niemeyer, Meerpohl, Stahl, Schratt, Glader, Backer, Wong, Nathan, Beggs and Sieff2006).
The identification of new splice variants in this study and the findings of their association with blastocyst formation provide a platform for further investigations of these genes. In addition, given that the splice variants of the investigated genes were differentially expressed in the embryo groups and that alternative splicing is a key factor in phenotypic complexity (Blencowe, Reference Blencowe2006), it is critical to generate a catalog of all splice variants of the bovine genome. Such a database would help us understand some of the mechanisms leading to phenotypic variation in cattle.
In conclusion, in this study we have investigated the association of novel transcripts in bovine embryos with early embryonic development. Some of the novel transcripts showed sequence similarity to human and cattle genes known to affect differentiation, growth and development. Moreover, we have characterized some alternatively spliced genes that showed differential alternative splicing between degenerate and normally developed embryos. One limitation in this study was the use of one RNA-Seq sample for each embryo population. However, for validation of RNA-Seq results, three biological replicate pools of blastocysts and degenerate embryos were used in qRT-PCR analysis. It is well recognized that despite the significant advances in IVF procedures in both humans and cattle, the success rates are still low. Indeed, a major barrier to improved IVF outcomes is the inability to reliably predict embryo competence for optimal embryo transfer. Furthermore, gene-specific expression profiles can provide a better understanding of the development of preimplantation embryos, which in turn would enhance fertility success. As such, genes reported in this study could be used as potential markers for embryo quality and development.
Acknowledgments
This study was supported by USDA Hatch grant WIS-142-PRJ17PH from the University of Wisconsin – Madison. Bao Zhang was supported by the China Scholarship Council.
Supplementary Material
The supplementary material referred to in this article is available online at http://www.journals.Cambridge.org/anm
