Inter-genus gene expression analysis in livestock fibroblasts using reference gene validation based upon a multi-species primer set

Abstract

Quantitative reverse transcription PCR (RT-qPCR) remains as an accurate approach for gene expression analysis but requires labor-intensive validation of reference genes using species-specific primers. To ease such demand, the aim was to design and test a multi-species primer set to validate reference genes for inter-genus RT-qPCR gene expression analysis. Primers were designed for ten housekeeping genes using transcript sequences of various livestock species. All ten gene transcripts were detected by RT-PCR in Bos taurus (cattle), Bubalus bubalis (buffaloes), Capra hircus (goats), and Ovis aries (sheep) cDNA. Primer efficiency was attained for eight reference genes using B. taurus—O. aries fibroblast cDNA (95.54–98.39%). The RT-qPCR data normalization was carried out for B. taurus vs. O. aries relative gene expression using Bestkeeper, GeNorm, Norm-finder, Delta CT method, and RefFinder algorithms. Validation of inter-genus RT-qPCR showed up-regulation of TLR4 and ZFX gene transcripts in B. taurus fibroblasts, irrespectively of normalization conditions (two, three, or four reference genes). In silico search in mammalian transcriptomes showed that the multi-species primer set is expected to amplify transcripts of at least two distinct loci in 114 species, and 79 species would be covered by six or more primers. Hence, a multi-species primer set allows for inter-genus gene expression analysis between O. aries and B. taurus fibroblasts and further reveals species-specific gene transcript abundance of key transcription factors.

About the Authors:

Marcelo T. Moura

Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Writing – original draft, Writing – review & editing

* E-mail: [email protected]

Affiliation: Departamento de Medicina Veterinária, Universidade Federal Rural de Pernambuco, Pernambuco, Brazil

ORCID logo http://orcid.org/0000-0003-3808-6607

Roberta L. O. Silva

Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Methodology, Project administration, Writing – original draft

Affiliation: Departamento de Genética, Universidade Federal de Pernambuco, Pernambuco, Brazil

Pábola S. Nascimento

Roles Formal analysis, Methodology, Resources

Affiliation: Departamento de Medicina Veterinária, Universidade Federal Rural de Pernambuco, Pernambuco, Brazil

José C. Ferreira-Silva

Roles Formal analysis, Investigation, Methodology, Resources

Affiliation: Departamento de Medicina Veterinária, Universidade Federal Rural de Pernambuco, Pernambuco, Brazil

Ludymila F. Cantanhêde

Roles Data curation, Formal analysis, Investigation, Methodology, Resources

Affiliation: Departamento de Medicina Veterinária, Universidade Federal Rural de Pernambuco, Pernambuco, Brazil

Ederson A. Kido

Roles Data curation, Formal analysis, Investigation, Methodology, Resources

Affiliation: Departamento de Genética, Universidade Federal de Pernambuco, Pernambuco, Brazil

Ana M. Benko-Iseppon

Roles Conceptualization, Funding acquisition, Investigation, Methodology, Project administration, Supervision, Writing – review & editing

Affiliation: Departamento de Genética, Universidade Federal de Pernambuco, Pernambuco, Brazil

Marcos A. L. Oliveira

Roles Funding acquisition, Project administration, Resources, Supervision, Writing – review & editing

Affiliation: Departamento de Medicina Veterinária, Universidade Federal Rural de Pernambuco, Pernambuco, Brazil

Introduction

Cellular states are the result of differential gene expression that arises from stepwise instructive factors during development [1,2]. Alternatively, disease states or phenotypic perturbation of cellular transcriptional programs by environmental or/and genetic factors can be readily identified by gene expression changes [3,4]. Appropriate assessment of transcript abundance is thus paramount for identification of such key cell state regulators and gene expression signatures [3,5].

Quantitative reverse transcription PCR (RT-qPCR) stands out as a method of choice for relative gene expression analysis [6–8]. It allows fast, accurate real-time measurement of gene expression, even for genes with low-transcript abundance [7]. Furthermore, RT-qPCR is also recurrently used to validate gene expression data stemming from other mRNA-based technologies, such as cDNA microarrays and RNA sequencing [6,9,10].

The RT-qPCR data can be presented as an absolute measurement (e.g. copy number per cell) or at a relative expression level basis between two or more experimental conditions [11]. Accurate measurements of relative gene expression by RT-qPCR require internal control genes (i.e., reference genes—RGs) to counterbalance possible bias from technical or experimental factors, such as RNA extraction methods, mRNA integrity, cell-type, developmental stage, ploidy level, and environmental conditions, among others [12–17].Thus adequate choice and validation of RGs remains as one of the most critical steps for the accuracy of the RT-qPCR assay and is mandatory for each experimental design [8,14,18–20].The algorithms built in software packages for RG selection play a major role in the stability value that is calculated for each candidate RG. Moreover, RG stability can be determined by its overall gene expression variation [21], intra-group expression variation [22], based on pair-wise CT analysis [23], by accounting for mRNA integrity [24], or combining multiple algorithms at once [25]. Another substantial progress was made by adopting an RG index rationality [21,26], since no single locus can be efficiently used as an internal control gene across multiple tissues or experimental conditions [27].

With the current relatively ease to generate large-scale gene expression datasets and the increasing demand for high-throughput or automated platforms [6,28,29], the development of tools for various species or technologies can become extremely advantageous. The use of a universal or multi-species primer set for PCR is an attractive approach to shortcut efforts to species-specific primer designing and testing [30–33]. To the best of our knowledge, despite the availability of a report for Capra hircus and Bos taurus [13], there is no multi-species primer set (MSPS) validated for intra-specific relative gene expression for several livestock species or for inter-genus RT-qPCR analysis. Therefore, the aim was to design a MSPS for selected livestock species that led to validation of RGs in Ovis aries and B.taurus primary cells and for inter-genus relative gene expression via RT-qPCR.

Material and methods

Ethics approval

The project was formally approved by the Ethics Commission (CEUA) at the Universidade Federal Rural de Pernambuco (UFRPE) under the license 031/2016, in full accordance with Brazilian law requirements for animal experimentation.

Primary somatic cell culture

Primary fibroblast cultures were established and propagated as previously described by Moura et al. [34]. Ear skin biopsies from artiodactyla males of B. taurus (cattle), Bubalus bubalis (buffaloes), Capra hircus (goats), and O. aries (sheep) were collected at local abattoirs (Pernambuco, Brazil), washed in saline solution (0.9% NaCl) with antibiotics at 35°C, and further washed in 70% ethanol, besides hair removal with a blade. Skin samples were transported to the laboratory in D-PBS containing amphotericin (Fungizone) and antibiotics at 4°C. Skin explants were obtained by slicing each ear biopsy and placing it on 35 or 60 mm Petri dishes. Explants were cultured in high glucose DMEM (Gibco) containing 10% FBS and antibiotics. Fibroblasts were passaged when dishes became confluent and in seven-day intervals afterwards. Early passage fibroblasts (passages 2 and 3) were used for total RNA extraction.

Ovaries of O. aries were collected at local abattoirs (Pernambuco, Brazil), transported in saline solution (0.9% NaCl) with antibiotics at 35°C, within three hours after slaughter. Cumulus-oocyte complexes (COC) were aspirated from 3–6 mm follicles with an 18G needle into 10 mL syringes, respectively. Follicular content was deposited in Petri dishes in HEPES buffered in vitro maturation medium (H-IVM), consisting of TCM-199 with Hank’s salts, 10% fetal bovine serum (FBS), 50 IU mL-1 heparin, and 50 μg mL-1 gentamicin sulfate. Cumulus cells were isolated from oocytes by gentle pipetting and used for total RNA extraction.

Total RNA extraction

Primary fibroblasts and cumulus cells (>1.0 x 106) of selected livestock species (B. taurus, B. bubalis, C. hircus, and O. aries) were washed in 1X PBS (1.15 g Na2HPO4, 0.2g KH2PO4, 8.0g NaCl, and 0.2g KCl), snap-frozen in liquid nitrogen, and stored at– 80°C. The total RNA extraction was carried out using Reliaprep RNA Cell Miniprep (Promega), as described by the manufacturer. The total RNA quality (260–280 and 260–230 ratios) and amount was determined by both Nanodrop 2000 C and Qubit (Thermo Scientific), respectively. The RNA samples were evaluated for their integrity by electrophoresis in 1.0% agarose gels, with 70 V, 120 A after a 40-min. run.

The cDNA synthesis

The reverse transcription (RT) reaction was carried out immediately after total RNA extraction. The procedure was performed with Quantitect RT kit (Qiagen), as described by the manufacturer. Possible residual genomic DNA was removed by the genomic DNA (gDNA) elimination reaction (2 μL gDNA wipeout buffer 7x, 1 μg total RNA and a final volume of 14 μL), incubated at 42°C for three min., and transferred to 4°C immediately. The gDNA elimination reaction was added to the RT reaction (1 μL Quantiscript RT, 4 μL Quantiscript RT buffer 5X, and 1 μL RT primer mix) and incubated at 42°C for 30 min. Moreover, the samples were incubated at 95°C for three min. and stored at—20°C.

Primer design

Primers were designed for the MSPS using Primer-BLAST [35] and further selected using Primer3plus [36] based on B. taurus, C. hircus, or O. aries EST sequences available in the GenBank (Table 1). Primers were designed for amplicons of 70–200 base pairs (bp), 40–60% GC content, 19–23 bases in size, to contain GC clamps while avoiding primer-dimer formation [36]. Only primers that were predicted to amplify gene transcripts from at least three livestock species (i.e., B. taurus, C. hircus, and O. aries) as assessed by Primer-BLAST were used in the experiment.

[Figure omitted. See PDF.]

Table 1. Multi-species primer set to validate reference genes for inter-genus quantitative reverse transcription PCR (RT-qPCR) in Ovis aries and Bos taurus fibroblasts.

https://doi.org/10.1371/journal.pone.0221170.t001

Reverse transcription polymerase chain reaction (Rt-PCR)

The RT-PCR reactions for the MSPS using livestock fibroblast cDNA (i.e., B. taurus, B. bubalis, C. hircus, and O. aries) were carried out in a 96 reaction thermal cycler (Techne), as previously described by Moura et al. [34]. The reaction was composed of 1 μL cDNA (non-diluted samples), 0.6 μL of each primer (2.5 μM), 2 μL 10 X PCR Buffer, 0.2 μL dNTP (10 mM), 0.6 μL 50 mM MgCl2, 0.1 μL Taq Polymerase (LGC Biotecnologia), and 14.9 μL ultra-pure H2O to a final volume of 20 μL. The cDNA was denatured, and Taq polymerase was activated at 94°C for 5 min. Moreover, 35 PCR cycles were carried out with the following conditions: denaturation at 94°C for 1 min., annealing at 58°C (C. hircus and O. aries) or 60°C (B. taurus and B. bubalis) for 1 min., and extension at 72°C for 1 min. The final extension step was performed at 72°C for 10 minutes. Amplicon visualization was carried out by electrophoresis in 0.5X TBE buffer and 1.5% agarose gels with 70 V, 120 A for 20 min.

Quantitative reverse transcription PCR (RT-qPCR)

The RT-qPCR reactions were performed in Line Gene 9660 FQD-96A (Bioer), using SYBR Green detection system. The reaction composed of 1 μL cDNA, 5 μL 2x Go Taq qPCR Master Mix (Promega), 0.6 μL primers (2.5 μM), and 3.4 μL ultra-pure H2O to a final volume of 10 μL. Both fibroblast and cumulus cell cDNA were used as diluted samples (1:10). The reactions were carried out under the following conditions: initial denaturation at 95°C for 2 min., 40 cycles at 95°C for 15 seconds, and 58°C for 60 seconds (O. aries cumulus cells cDNA) or 60°C for 60 seconds (O. aries fibroblast cDNA or B. taurus fibroblast cDNA) for relative gene expression analysis. Samples of bulk fibroblast cDNA (1:1 O. aries—B. taurus) were initially used for determining primer efficiency in a similar fashion to as previously described for sheep-goat hybrid embryos [37], followed by species-specific primer efficiency testing (using non-mixed cDNA) and relative gene expression analysis (S1 Fig). Melting curves were analyzed from 65 to 95°C for 20 min. after the initial 40 cycles. Expression levels of candidate genes were evaluated by the cycle of quantification (Cq) during the exponential phase (log) of the PCR. Primer amplification efficiency (E = 10−1 / slope), correlation coefficient (R2), and interception (y) were determined by the standard curve method using cDNA serial dilutions: 100 (non-diluted samples), 10−1, 10−2, 10−3, and 10−4 [38]. The Minimum Information for Publication of Quantitative Real-Time PCR Experiments (MIQE) guidelines was followed to increase both the transparency and reliability of the results generated in the study [12,14,19,20].

The relative gene expression was obtained with three biological replicates and three technical replicates. The relative gene expression levels were evaluated with the REST tool (version 2.0.13) that relies on the 2-ΔΔCT method [11,39], which is based on pairwise comparisons using randomization and bootstrapping techniques—Pairwise Fixed Reallocation Randomization Test [40].

RT-qPCR data normalization

The BestKeeper algorithm uses raw Cq values as input [21], and calculates the standard deviation (SD) from Cq values for each RG. The RGs with >1.0 SD are considered inconsistent and should be removed from the analysis. By a Pearson correlation coefficient (geometric mean of Cq values of candidate genes) the BestKeeper (version 1) identifies the most stable RGs, defining the RG index for RT-qPCR data normalization [21].

The model-based approach to estimation of expression variation by the NormFinder algorithm ranks the most suitable genes from the candidate RG set [22] while considering both intragroup and intergroup variations. Further, the NormFinder (version 0.953) identifies the most stable RG and the two most appropriate RG by pair-wise comparisons. In turn, the GeNorm algorithm (version 3.5) calculates the RG stability measurement (M) as the average pairwise variation of each RG with all the other RGs [26]. It also enables the elimination of less stable RGs and reassessment of M-values, resulting in the final ranking of RGs. The RGs with lower M-value hold more stable expression among samples [26]. Moreover, GeNorm algorithm also calculates the RG index after inclusion of additional RGs to the most stable RG pair initially found (Pairwise variation—V).

The delta-CT method was also used for ranking candidate RGs [23]. This algorithm is based on a simple pairwise delta-CT measurement in each sample. The readout of the delta-CT method is a RG ranking from the most stable to the least stable candidate RG. The raw RT-qPCR data was also subject to the RefFinder algorithm to combine the results from the aforementioned software packages [25], which is available online at the following link (https://www.heartcure.com.au/for-researchers/).

In silico analysis of primer specificity and candidate reference gene expression in somatic tissues

RNA sequencing (RNA-Seq) data for each candidate RG of both B. taurus and O. aries was retrieved from gene expression atlas at the European Bioinformatics Institute (https://www.ebi.ac.uk/gxa/home). Only RNA-Seq experiments that contained common tissues in both species were analyzed and described in a side-by-side fashion.

In silico analysis of primer specificity was carried out using Primer-BLAST [35]. Each primer of the MSPS was subject to BLAST to mammalian transcriptomes in the GenBank database (Tax-id: 40674) or selected non-mammalian species (Anolis carolinensis, Danio rerio, Drosophila melanogaster, Gallus gallus, and Xenopus laevis). A maximum of two mismatches in each primer was tolerated for positive hits. The Primer-BLAST search was performed in 05.24.2018. All results were curated manually.

Results

Designing a multi-species primer set and testing its specificity in selected livestock species by RT-PCR and RT-qPCR

Total RNA was initially retrieved for two-step RT-qPCR (S2 Fig). The experimental design for primer design and validation (S3 Fig) was effective at obtaining primers for all ten candidate gene transcripts (Table 1). Based on in silico analysis, all primers expected to lead to specific PCR reactions in B. taurus, C. hircus,and O. aries. The primer design also allowed primers annealing at conserved sequences (S3 Fig), covering multiple deposited sequences of each loci in B. taurus, B. bubalis, C. hircus, and O. aries (S3 Fig).

All ten gene transcripts of the MSPS detected in O. aries cDNA from primary fibroblasts, displaying single amplicons with their expected sizes (Fig 1). Under identical PCR conditions, the MSPS was also detected in C. hircus fibroblasts samples (Fig 1). Moreover, adjusting annealing temperature allowed the detection of single amplicons for all ten gene transcripts of the MSPS in B. taurus and B. bubalis fibroblast cDNA (Fig 1). Additionally, the MSPS was also detected as absolute RT-qPCR reactions (B. taurus, B. bubalis, C. hircus, and O. aries), as demonstrated by the merged melting curves from qPCR reactions of these four species. Of note, subtle differences in melting curves may have arisen from varying isoform amplification or species-specific amplification efficiency (S4 Fig).

[Figure omitted. See PDF.]

Fig 1. Detection of gene transcripts from the multi-species primer set (MSPS) in C. hircus, O. aries, B. taurus, and B. Bubalis fibroblast cDNA by reverse transcriptase polymerase chain reaction (RT-PCR).

Molecular weight: 100 base pair (bp).

https://doi.org/10.1371/journal.pone.0221170.g001

Primer efficiency for RT-qPCR using a multi-species primer set and cDNA from O. aries and B. taurus primary somatic cells

Primer efficiency was determined using B. taurus—O. aries fibroblast cDNA (Table 2). From the ten primers, eight were efficient (95.54% - 98.39 while SDHA and TBP did reach this threshold, and were therefore removed from the data representation (Table 2), being not further investigated. The correlation coefficient varied from—0.990 to—0.999, with slopes from—3.37 to—3.43, and Y intercepts from 26.74 to 31.33 (Table 2). To exclude the potential of bulked cDNA-based bias, primer efficiency was also tested for a subset of primers using non-mixed O. aries and B. taurus fibroblast cDNA (S4 Table). Furthermore, primer efficiency in O. aries varied from 90.69% - 103.07%, while B. taurus primer efficiency ranged from 92.26–102.04%. The variation for correlation coefficient (- 0.990 to—0.999), slope (- 3.25 to—3.57), and Y intercepts (24.74 to 32.59) were slightly higher using non-mixed cDNA in comparison to the bulk cDNA approach (S4 Table). In turn, the determination of primer efficiency using O. aries cumulus cDNA found that nine primers were efficient (90.69% - 110.76%), with the exception of H3F3A (S1 Table). Also using sheep cumulus cDNA, the correlation coefficient varied from—0.995 to—0.999, with slopes from—3.09 to—3.57, and Y intercepts from 24.74 to 31.83. All primers of the MSPS were efficient in at least one experimental condition or cell type, thus suggesting that their efficiency was dependent on cell type-specific unidentified factors rather than limitations stemming in primer design.

[Figure omitted. See PDF.]

Table 2. Primer efficiency, coefficient correlation, slope, and Y intercept derived from the standard curve of each candidate reference gene from a multi-species primer set using O. aries—B. taurus fibroblast cDNA via RT-qPCR assay.

https://doi.org/10.1371/journal.pone.0221170.t002

RT-qPCR data normalization

The BestKeeper algorithm ranked the seven RG candidates by their Pearson coefficient (Table 3). However, both RPL19 and YWHAZ were considered inconsistent due to their high standard deviation (1.18 and 1.16) and removed from the analysis (Table 3). Based on the Pearson coefficient, ATP1A1 and UBB showed higher values than GAPDH, PPIA, and ACT, thus leading to a more stable RG index as a pair.

[Figure omitted. See PDF.]

Table 3. Standard deviations (SD) and Pearson’s correlations obtained from RT-qPCR data using the BestKeeper algorithm.

https://doi.org/10.1371/journal.pone.0221170.t003

The RT-qPCR data normalization using NormFinder algorithm found ATP1A1 as the most stable RG (0.011). The best RG pair using the NormFinder was the same as the BestKeeper algorithm, namely ATP1A1 and UBB (0.010). Moreover, the RG ranking by NormFinder was substantially different from BestKeeper (Table 4), although three out of the four most stable were shared by both algorithms.

[Figure omitted. See PDF.]

Table 4. Stability values derived from RT-qPCR data using B. taurus and O. aries fibroblasts determined by NormFinder and GeNorm algorithms.

https://doi.org/10.1371/journal.pone.0221170.t004

The GeNorm algorithm ranked the RGs in a similar way as the NormFinder (Table 4). After the GeNorm RG pairwise analysis, RPL19 and YWHAZ was the most stable RG pair (Fig 2). The pairwise variation (V) analysis showed that increasing the number of RG to three and four, also increases the RT-qPCR data normalization efficiency (Fig 2).

[Figure omitted. See PDF.]

Fig 2.

Average gene expression stability values (M) of seven potential reference genes from the multi-species primer set (MSPS) for normalization of inter-genus RT-qPCR data based on the GeNorm analysis (A). Pairwise variation—V and inter-genus analysis of seven potential reference genes for RT-qPCR data normalization (B).

https://doi.org/10.1371/journal.pone.0221170.g002

The delta-CT method ranked ATP1A1, RPL19, YWHAZ, and UBB as the four best ranked RGs (Fig 3). It is noteworthy that RG ranking by the delta-CT method was substantially different from GeNorm and NormFinder algorithms (Fig 3), but more in agreement with BestKeeper. Both ACT and GAPDH were the two least stable candidate RGs (Fig 3).

[Figure omitted. See PDF.]

Fig 3. Reference gene ranking by the delta-CT method for normalization of B. taurus–O. aries inter-genus RT-qPCR.

https://doi.org/10.1371/journal.pone.0221170.g003

The RefFinder algorithm ranked ATP1A1, RPL19, UBB, and YWHAZ as the top four most stable RGs (Fig 4), thus very closely overlapping the delta-CT method (Fig 3). The exceptions to this rule are YWHAZ and UBB, which ranked interchanging at the third and fourth between these two algorithms.

[Figure omitted. See PDF.]

Fig 4. Reference gene ranking by the RefFinder algorithm for normalization of B. taurus–O. aries inter-genus RT-qPCR.

https://doi.org/10.1371/journal.pone.0221170.g004

Validation of inter-genus RT-qPCR using Toll-like receptor 4 (TLR4) and Zinc finger protein, X-linked (ZFX) gene transcripts

Validation of inter-genus relative gene expression via RT-qPCR was carried out for TLR4 and ZFX using four RGs (ATP1A1, RPL19, UBB, and YWHAZ), under varying normalization regimens. Using the two most stable RGs found by BestKeeper and NormFinder algorithms (i.e., ATP1A1 and UBB), ZFX was found to be induced in B. taurus by 2.98 fold, while TLR4 was also up-regulated by 7.59 fold (Fig 5). The RT-qPCR normalization using the best RG pair found by the GeNorm algorithm (ATP1A1 and RPL19), showed that candidate genes display a milder but consistently higher expression in B. taurus fibroblasts (ZFX: 1.43 fold and TLR4: 2.63 fold). The addition of a third and fourth RGs to the best RG pair found by GeNorm increased the relative gene expression difference between species (Fig 5). Collectively, the results suggest that both candidate genes are up-regulated in B. taurus fibroblasts in comparison to their O. aries counterparts.

[Figure omitted. See PDF.]

Fig 5.

Relative expression of Toll-like receptor 4—TLR4 (A) and Zinc finger protein, X-linked—ZFX (B) gene transcripts determined by the REST software using the 2-ΔΔCT methodology in cDNAs from O. aries and B. taurus fibroblasts via inter-genus RT-qPCR analysis. A: ATPase Na+/K+ transporting subunit alpha 1 (ATP1A1). R: Ribosomal protein L19 (RPL19). U: Ubiquitin B (UBB). Y: Tyrosine 3-monooxygenase/tryptophan 5-monooxygenase activation protein zeta (YWHAZ). Blue line: NormFinder and Bestkeeper-based RT-qPCR normalization. Green line: GeNorm RT-qPCR normalization using 2, 3, and 4 genes. Purple line: Other combinations of reference genes for RT-qPCR normalization.

https://doi.org/10.1371/journal.pone.0221170.g005

In silico approach for assessing the multi-species primer set for several B. taurus and O. aries tissues and for other mammals

To better estimate the potential of the MSPS for normalizing RT-qPCR data in B. taurus and O. aries, RNA-Seq data was retrieved for seven organs (tissues) from an available gene expression atlas (S5 Fig). The most stable RGs above (ATP1A1, RPL19, UBB, and YWHAZ) are detected in most tissues, except muscle of B. taurus (S5 Fig).

A broader primer-BLAST search of each primer from the MSPS for mammalian species showed that they would amplify at least two gene products in 114 species, from 24 orders and 96 genera (S2 Table). A clear distinction is made from the number of species covered by a varying number of primers (Fig 6). Most species are covered by five to seven primers, followed by species attained by eight to ten primers (Fig 6). Only eight species are covered by all ten primers, but the number increases to 79 species for six or more primers (Fig 6).

[Figure omitted. See PDF.]

Fig 6. Multi-species primer set (MSPS) specificity in mammals.

Number of species covered by varying primers, MSPS coverage by the number of primers, and primer-specific specificity.

https://doi.org/10.1371/journal.pone.0221170.g006

Based on the Primer-BLAST search, the expected rate of positive specific reactions was 98.16% (749 of 763 potential reactions), with 1.57% (12 of 763) of non-specific reactions (i.e., gene transcript of interest plus other unintended transcripts), and 0.26% unspecific reactions (2 of 763 potential reactions). The in silico primer specificity among the MSPS in all retrieved mammalian species is expected to vary from 93.69% of H3F3A to 100.00% for several primers (ATP1A1, SDHA, TBP, UBB, and YWHAZ).

Most species found in the primer-BLAST search were representative of five orders (Artiodactyla, Rodentia, Primata, Carnivorae, and Chiroptera; S2 Table), thus representing 73.68% of the mammals. There were clear distinctions of primer specificity among orders and species (S6 Fig). As expected, most Artiodactyla species (ruminants) were covered by all primers (S2 Table), but with some species-specific exceptions, mostly for UBB (S6 Fig). Curiously, Bos indicus displayed polymorphisms that are predicted to not allow SDHA and TBP gene product amplifications, in contrast to B. taurus and Bos mutus (S2 Table). Several non-ruminant artiodactyla species were also covered by the MSPS as determined by primer-BLAST (S2 Table). Primate transcriptomes were mostly covered by six primers (ACTIN, ATP1A1, H3F3A, PPIA, RPL19, and YWHAZ), while rodents were generally covered by four (ACTIN, ATP1A1, H3F3A, and TBP), carnivores by eight (ACTIN, ATP1A1, GAPDH, H3F3A, RPL19, SDHA, TBP, and YWHAZ), and Chiroptera by four primers (ACTIN, GAPDH, PPIA, and RPL19). Finally, the assessment by primer-BLAST of other five non-mammalian species showed that one to three primers are expected to amplify gene transcripts (S3 Table).

Discussion

Gene expression analysis has been a fruitful tool for understanding the establishment and maintenance of cellular states [1–2]. However, less attention has been drawn to RT-qPCR conditions, further leading to data of lower reproducibility [8,14,19]. Within the requirements of reliable RT-qPCR, RG selection remains paramount and must be carried out for each experimental setup [13–18].

The strategy of designing primers at conversed sequences for selected livestock species allowed the development of a MSPS for potential use in the RT-qPCR assay. Although “universal” primers have been described in different contexts [30–33], to the best of our knowledge, a MSPS remained to be described for RT-qPCR applicable to multiple order-specific species and for inter-genus gene expression analysis. There are several evidences that the here described MSPS is useful for the livestock species; firstly, all ten primers detected as single amplicons by conventional RT-PCR in four different genera. Secondly, all ten primers efficient in O. aries—B. taurus bulk cDNA, O.aries and B. taurus fibroblast cDNA or O. aries cumulus cDNA, suggesting that primer efficiency was dependent upon cell-type specific gene transcript abundance. Furthermore, the multi-species mixed cDNA approach has been validated under other experimental contexts [37,41]. Thirdly, RT-qPCR normalization was attained in B. taurus—O. aries inter-genus RT-qPCR as assessed by five algorithms. Fourthly, some primers have been used for intra-specific RT-PCR in eggs and early embryos of livestock species [42,43].

Different RG sets were found for O. aries—B. taurus inter-genus RT-qPCR-based relative gene expression normalization using three algorithms. The GeNorm normalization using RPL19 and YWHAZ showed a substantial difference to the NormFinder ATP1A1 and UBB-pairwise index, most likely due to methodological differences. However, the addition of the third and fourth RGs to the GeNorm RG pair demonstrated a parallel between normalization software packages. Both methods (delta-CT and RefFinder) generated similar RG rankings but did not deviate far from the other aforementioned algorithms. Numerous strategies have been initially described for identifying a single RG [18], and later the RG index approach was introduced for more accurate RT-qPCR data normalization [21,26]. Similar reports for RG validation used various numbers of candidate gene transcripts and algorithms [13–17,44]. More importantly, some reports have found similar RG sets using different software packages [45], while other have reported contrasting RG sets [46], particularly using the most widely used algorithms [21,26,44].

The widely used ACT and GAPDH genes proved to be within the most unstable ones, as it is becoming evident in many other experimental conditions and organisms [27,47], further highlighting the importance of careful and unbiased RG selection [21–26]. Novel not yet tested candidate RG have also proven to be more stable than the widely used ones [47], at least under some circumstances. The critical selection of candidate RG also has the capability to diminish or avoid potential selection of co-regulated RGs, which may weaken RT-qPCR data normalization [48]. Primer efficiency variation is another critical experimental parameter that may affect RT-qPCR data normalization [48], as observed for the gene of choice TLR4, as described above Nonetheless, no primer efficiency variation was found, since RG primer efficiency was attained using single-species or bulked cDNA of B. taurus—O. aries from fibroblasts, as previously described elsewhere [37,38]. Despite the fact that algorithms indicated different RG indexes or rankings, it became clear that all results from such software packages converged to similar results. Collectively, these findings demonstrate that different RG indexes can be successfully used for inter-genus RT-qPCR data normalization. Further evidence suggests that the retrieved RG described above may indeed be stable across other tissues, as demonstrated by RNA-Seq data from seven somatic tissues.

Both TRL4 and ZFX candidate genes were up-regulated in B. taurus fibroblasts in comparison to their O. aries counterparts. Toll-like receptors (TLRs) are transmembrane-spanning proteins and are part of the innate immune response, displaying high species-specific, tissue-specific, and individual variation in transcript availability [49,50], thus reinforcing the results described above. The ZFX is an X-linked gene that encodes a transcription factor expressed throughout development [51–53]. Since ZFX is a highly conserved protein across mammals [54], differences in transcript availability between O. aries and B. taurus are possibly caused by differences in cis-regulatory sequences or varying trans-activating proteins. Gene expression profiling across species is an instructive tool for understanding species-specific gene expression signatures [3,5]. However, less attention has been drawn to the quantitative measurements of gene expression among species [13,16], particularly in a direct, unbiased fashion [55,56]. Collectively, these results suggest that fibroblasts may holdspecies-specific gene stoichiometry as described in other systems [57,58], and thus may have biological significance as for cell potency, lineage specification, and cellular reprogramming [57–59]. However, further evidence is required for full appreciation of this variation at transcriptomic level across species. Alternatively, the MSPS described here could also be used for other applications such as diagnostics [60], genetic characterization [61], nucleus-cytoplasm interactions [62], and interspecies cell-fusion or nuclear transfer experiments [56,63].

An in silico search in mammals using primer-BLAST showed that MSPS is predicted to amplify at least two gene transcripts of over a hundred species. More importantly, a total of around eighty species would be covered by six or more primers from the MSPS. These results are striking by the fact that 24 mammalian orders were represented in the initial in silico screening; with around half of such species are rodents or primates, including many common animal models for developmental biology and biomedical research (e.g., rat, rhesus, and gorilla). The MSPS can also be improved for order-specific RT-qPCR by the design of a small subset of primers. Additionally, the results reinforce that the potential of the MSPS has been underestimated and could be potentially extended to other classes of gene transcripts for cross-species or inter-genus gene expression profiling via RT-qPCR. In sum, the design of a multi-species primer set allows for inter-genus gene expression analysis between O. aries and B. taurus primary cells and further reveals species-specific gene transcript abundance of key transcription factors.

Supporting information

[Figure omitted. See PDF.]

S1 Table. Primer efficiency, coefficient correlation, slope, and Y-intercept derived from the standard curve of each candidate reference gene from the multi-species primer set (MSPS) using O. aries cumulus cDNA via RT-qPCR assay.

https://doi.org/10.1371/journal.pone.0221170.s001

(PDF)

S2 Table. In silico analysis of multi-species primer set (MSPS) specificity in mammalian orders determined by primer-BLAST.

https://doi.org/10.1371/journal.pone.0221170.s002

(PDF)

S3 Table. In silico analysis of multi-species primer set (MSPS) specificity in selected non-mammalian species determined by primer-BLAST.

https://doi.org/10.1371/journal.pone.0221170.s003

(PDF)

S4 Table. Primer efficiency, coefficient correlation, slope, and Y-intercept derived from the standard curve of each candidate reference gene from the multi-species primer set (MSPS) using non-mixed O. aries and B. taurus fibroblast cDNA via RT-qPCR assay.

https://doi.org/10.1371/journal.pone.0221170.s004

(DOCX)

S1 Fig. Description of cDNA type used for primer efficiency testing and relative gene expression analysis.

https://doi.org/10.1371/journal.pone.0221170.s005

(TIFF)

S2 Fig. Total RNA extraction efficiency and RNA 280/260 and 260/230 ratios from somatic cell samples of selected livestock species (B. taurus, B. bubalis, C. hircus, and O. aries) investigated in the experiment.

https://doi.org/10.1371/journal.pone.0221170.s006

(TIFF)

S3 Fig.

Experimental design (A), site of primer annealing in gene transcripts (B), and number of gene transcripts covered by the multi-species primer set (MSPS) in selected livestock species (B. taurus, B. bubalis, C. hircus, and O. aries).

https://doi.org/10.1371/journal.pone.0221170.s007

(TIFF)

S4 Fig. Detection of gene transcripts (melting curves) from the multi-species primer set (MSPS) in selected livestock species by quantitative reverse transcription PCR (RT-qPCR).

https://doi.org/10.1371/journal.pone.0221170.s008

(PDF)

S5 Fig. Relative expression of candidate reference genes in seven tissues of B. taurus and O. aries as determined by the gene expression atlas at the European Bioinformatics Institute website (https://www.ebi.ac.uk/gxa/home).

TPM: transcripts per million.

https://doi.org/10.1371/journal.pone.0221170.s009

(TIFF)

S6 Fig. Order-specific multi-species primer set (MSPS) specificity in mammals.

Number of species covered by varying primers, MSPS coverage by the number of primers, and primer-specific specificity.

https://doi.org/10.1371/journal.pone.0221170.s010

(TIFF)

Acknowledgments

We would like to acknowledge Sérgio Nascimento for assistance with cell culture. We also thank CNPq (National Council for Scientific and Technological Development) and CAPES (Coordination of Higher Level Personnel) for financial support and fellowships.

Citation: Moura MT, Silva RLO, Nascimento PS, Ferreira-Silva JC, Cantanhêde LF, Kido EA, et al. (2019) Inter-genus gene expression analysis in livestock fibroblasts using reference gene validation based upon a multi-species primer set. PLoS ONE 14(8): e0221170. https://doi.org/10.1371/journal.pone.0221170

References

1. Armond JW, Saha K, Rana AA, Oates CJ, Jaenisch R, Nicodemi M, et al. Stochastic model dissects cell states in biological transition processes. Sci Rep. 2014;4: 3692. pmid:24435049

2. Gonzales KA, Ng HH. Biological Networks Governing the Acquisition, Maintenance, and Dissolution of Pluripotency: Insights from Functional Genomics Approaches. Cold Spring HarbSymp Quant Biol. 2015;80: 189–198.

3. Montgomery SB, Dermitzakis ET. The resolution of the genetics of gene expression. Hum Mol Genet. 2009;18: R211–R215. pmid:19808798

4. Burga A, Lehner B. Predicting phenotypic variation from genotypes, phenotypes and a combination of the two. Curr Opin Biotechnol. 2013;24: 803–809. pmid:23540420

5. Wittkopp PJ. Variable gene expression in eukaryotes: a network perspective. J Exp Biol. 2007;210: 1567–1575. pmid:17449821

6. VanGuilder HD, Vrana KE, Freeman WM. Twenty-five years of quantitative PCR for gene expression analysis. Biotechniques. 2008;44: 619–626. pmid:18474036

7. Derveaux S, Vandesompele J, Hellemans J. How to do successful gene expression analysis using real-time PCR. Methods. 2010;50: 227–230. pmid:19969088

8. Bustin S, Nolan T. Talking the talk, but not walking the walk: RT-qPCR as a paradigm for the lack of reproducibility in molecular research. Eur J Clin Invest. 2017;47: 756–774. pmid:28796277

9. Everaert C, Luypaert M, Maag JLV, Cheng QX, Dinger ME, Hellemans J, et al. Benchmarking of RNA-sequencing analysis workflows using whole-transcriptome RT-qPCR expression data. Sci Rep. 2017;7: 1559. pmid:28484260

10. Melo EO, Cordeiro DM, Pellegrino R, Wei Z, Daye ZJ, Nishimura RC, et al. Identification of molecular markers for oocyte competence in bovine cumulus cells. Anim. Genet. 2017;48: 19–29. pmid:27650317

11. Schmittgen TD, Livak KJ. Analyzing real-time PCR data by the comparative C T method. Nat Protoc. 2008;3: 1101–1108. pmid:18546601

12. Bustin SA.Why the need for qPCR publication guidelines?—The case for MIQE, Methods. 2010;50: 217–226. pmid:20025972

13. Bonnet M, Bernard L, Bes S, Leroux C.Selection of reference genes for quantitative real-time PCR normalisation in adipose tissue, muscle, liver and mammary gland from ruminants. Animal. 2013;7: 1344–1353. pmid:23552195

14. Edmunds RC, McIntyre JK, Luckenbach JA, Baldwin DH, Incardona JP. Toward enhanced MIQE compliance: reference residual normalization of qPCR gene expression data. J Biomol Tech. 2014;25: 54–60. pmid:24982597

15. Gu Y, Shen X, Zhou D, Wang Z, Zhang N, Shan Z, et al. Selection and Expression Profiles of Reference Genes in Mouse Preimplantation Embryos of Different Ploidies at Various Developmental Stages. PLoS One. 2014;9: e98956. pmid:24927500

16. Macabelli CH, Ferreira RM, Gimenes LU, De Carvalho NA, Soares JG, Ayres H, et al. Reference gene selection for gene expression analysis of oocytes collected from dairy cattle and buffaloes during winter and summer. PLoS One. 2014;9: e93287. pmid:24676354

17. Wang YK, Li X, Song ZQ, Yang CX. Methods of RNA preparation affect mRNA abundance quantification of reference genes in pig maturing oocytes. Reprod Domest Anim. 2017;52: 722–730. pmid:28407308

18. Huggett J, Dheda K, Bustin S, Zumla A. Real-time RT-PCR normalisation; strategies and considerations. Genes Immun. 2005;6: 279–284. pmid:15815687

19. Taylor S, Wakem M, Dijkman G, Alsarraj M, Nguyen MA. A practical approach to RT-qPCR-Publishing data that conform to the MIQE guidelines. Methods. 2010;50: S1–S5. pmid:20215014

20. Taylor SC, Mrkusich EM. The state of RT-quantitative PCR: firsthand observations of implementation of minimum information for the publication of quantitative real-time PCR experiments (MIQE). J Mol Microbiol Biotechnol. 2014;24: 46–52. pmid:24296827

21. Pfaffl MW, Tichopad A, Prgomet C, Neuvians TP. Determination of stable housekeeping genes, differentially regulated target genes and sample integrity: BestKeeper—Excel-based tool using pair-wise correlations. Biotechnol Lett. 2004;26: 509–515. pmid:15127793

22. Andersen CL, Jensen JL, Ørntoft TF. Normalization of real-time quantitative reverse transcription-PCR data: a model-based variance estimation approach to identify genes suited for normalization, applied to bladder and colon cancer data sets. Cancer Res. 2004;64: 5245–5250. pmid:15289330

23. Silver N, Best S, Jiang J, Thein SL. Selection of housekeeping genes for gene expression studies in human reticulocytes using real-time PCR. BMC Mol Biol. 2006;7: 33. pmid:17026756

24. Fleige S, Walf V, Huch S, Prgomet C, Sehm J, Pfaffl MW. Comparison of relative mRNA quantification models and the impact of RNA integrity in quantitative real-time RT-PCR. Biotechnol Lett. 2006;28: 1601–1613. pmid:16900335

25. Xie F, Xiao P, Chen D, Xu L, Zhang B. miRDeepFinder: a miRNA analysis tool for deep sequencing of plant small RNAs. Plant Mol Biol. 2012;80: 75–84.

26. Vandesompele J, De Preter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, et al. Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes. Genome Biol. 2002;3: research0034. pmid:12184808

27. Kozera B, Rapacz M. Reference genes in real-time PCR. J Appl Genet. 2013;54: 391–406. pmid:24078518

28. Gerrard G, Mudge K, Foskett P, Stevens D, Alikian M, White HE, et al. Fast-mode duplex qPCR for BCR-ABL1 molecular monitoring: innovation, automation, and harmonization. Am J Hematol. 2012;87: 717–720. pmid:22566190

29. Tan GW, Khoo AS, Tan LP. Evaluation of extraction kits and RT-qPCR systems adapted to high-throughput platform for circulating miRNAs. Sci Rep. 2015;5: 9430. pmid:25800946

30. Meiyu F, Huosheng C, Cuihua C, Xiaodong T, Lianhua J, Yifei P, et al. Detection of flaviviruses by reverse transcriptase-polymerase chain reaction with the universal primer set. Microbiol Immunol. 1997;41: 209–213. pmid:9130232

31. Pineau P, Henry M, Suspène R, Marchio A, Dettai A, Debruyne R, et al. A universal primer set for PCR amplification of nuclear histone H4 genes from all animal species. MolBiol Evol. 2005;22: 582–588.

32. Puech C, Dedieu L, Chantal I, Rodrigues V. Design and evaluation of a unique SYBR Green real-time RT-PCR assay for quantification of five major cytokines in cattle, sheep and goats. BMC Vet Res. 2015;11: 65. pmid:25889787

33. Jean A, Tardy F, Allatif O, Grosjean I, Blanquier B, Gerlier D. Assessing mycoplasma contamination of cell cultures by qPCR using a set of universal primer pairs targeting a 1.5 kb fragment of 16S rRNA genes. PLoS One. 2017;12: e0172358. pmid:28225826

34. Moura MT, Sousa RV, Leme LO, Rumpf R. Analysis of actinomycin D treated cattle oocytes and their use for somatic cell nuclear transfer. Anim Reprod Sci. 2008;109: 40–49. pmid:18162337

35. Ye J, Coulouris G, Zaretskaya I, Cutcutache I, Rozen S, Madden TL. Primer-BLAST: a tool to design target-specific primers for polymerase chain reaction. BMC Bioinformatics. 2012;13: 134. pmid:22708584

36. Hung JH, Weng Z. Designing Polymerase Chain Reaction Primers Using Primer3Plus. Cold Spring Harb Protoc. 2016;2016: pdb.prot093096.

37. Ma LB, Yang L, Zhang Y, Cao JW, Hua S, Li JX. Quantitative analysis of mitochondrial RNA in goat-sheep cloned embryos. Mol Reprod Dev. 2008;75: 33–39. pmid:17570506

38. Silva LR, Silva MD, Ferreira Neto JRC, de Nardi CH, Chabregas SM, Burnquist WL, et al. Validation of novel reference genes for reverse transcription quantitative real-time PCR in drought-stressed sugarcane. ScientificWorldJournal. 2014;2014: 357052. pmid:24987730

39. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2-ΔΔCT method. Methods. 2001;25: 402–408. pmid:11846609

40. Pfaffl MW, Horgan GW, Dempfle L. Relative expression software tool (REST) for groupwise comparison and statistical analysis of relative expression results in real-time PCR. Nucleic Acids Res. 2002;30: e36. pmid:11972351

41. Kleine M, Schrem H, Borlak J, Klempnauer J. Clinical versatility of porcine hepatocytes in the light of interspecies differences in cytochrome P450 regulation and expression.Xenotransplantation. 2008;15: 208–217. pmid:18957043

42. Moura MT, Silva RLO, Cantanhêde LF, Silva JB, Ferreira-Silva JC, Silva PG et al. Activity of non-canonical pluripotency-associated transcription factors in goat cumulus-oocyte complexes. Livest Sci. 2018;212: 52–56.

43. Silva PGC, Moura MT, Silva RLO, Nascimento PS, Silva JB, Ferreira-Silva JC, et al. Temporal expression of pluripotency-associated genes in sheep and cattle preimplantation embryos. Zygote. 2018;26: 270–278. pmid:30033902

44. Jacob F, Guertler R, Naim S, Nixdorf S, Fedier A, Hacker NF, Heinzelmann-Schwarz V. Careful selection of reference genes is required for reliable performance of RT-qPCR in human normal and cancer cell lines. PLoS One. 2013;8: e59180. pmid:23554992

45. Zhu W, Lin Y, Liao H, Wang Y. Selection of reference genes for gene expression studies related to intramuscular fat deposition in Capra hircus skeletal muscle. PLoS One. 2015;10: e0121280. pmid:25794179

46. Martínez-Giner M, Noguera JL, Balcells I, Fernández-Rodríguez A, Pena RN. Selection of internal control genes for real-time quantitative PCR in ovary and uterus of sows across pregnancy. PLoS One. 2013;8: e66023. pmid:23785467

47. Park SJ, Kim YH, Huh JW, Lee SR, Kim SH, Kim SU, et al. Selection of new appropriate reference genes for RT-qPCR analysis via transcriptome sequencing of cynomolgus monkeys (Macacafascicularis). PLoS One. 2013;8: e60758. pmid:23613744

48. Kaur R, Sodhi M, Sharma A, Sharma VL, Verma P, Swami SK, et al. Selection of suitable reference genes for normalization of quantitative RT-PCR (RT-qPCR) expression data across twelve tissues of riverine buffaloes (Bubalus bubalis). PLoS One. 2018;13: e0191558. pmid:29509770

49. Menzies M, Ingham A. Identification and expression of Toll-like receptors 1–10 in selected bovine and ovine tissues. Vet Immun Immunopathol. 2006;109: 23–30.

50. Dhanasekaran S, Kannan TA, Raj GD, Tirumurugaan KG, Raja A, Kumanan K. Differential expression of toll-like receptor mRNA in corneal epithelium of ruminants. Vet Ophthalmol. 2010;13: 270–274. pmid:20618807

51. Luoh SW, Bain PA, Polakiewicz RD, Goodheart ML, Gardner H, Jaenisch R, et al. ZFX mutation results in small animal size and reduced germ cell number in male and female mice. Development. 1997;124: 2275–2284. pmid:9187153

52. Galan-Caridad JM, Harel S, Arenzana TL,Hou ZE, Doetsch FK, Mirny LA, et al. Zfx controls the self-renewal of embryonic and hematopoietic stem cells. Cell. 2007;129: 345–357. pmid:17448993

53. Harel S, Tu EY, Weisberg S, Esquilin M, Chambers SM, Liu B, et al. ZFX controls the self-renewal of human embryonic stem cells. PLoS One. 2012;7: e42302. pmid:22879936

54. Poloumienko A. Cloning and comparative analysis of the bovine, porcine, and equine sex chromosome genes ZFX and ZFY. Genome. 2004;47: 74–83. pmid:15060604

55. Pei DS, Sun YH, Chen SP, Wang YP, Hu W, Zhu ZY. Zebrafish GAPDH can be used as a reference gene for expression analysis in cross-subfamily cloned embryos. Anal Biochem. 2007a;363: 291–293.

56. Pei DS, Sun YH, Chen SP, Wang YP, Hu W, Zhu ZY. Identification of differentially expressed genes from the cross-subfamily cloned embryos derived from zebrafish nuclei and rare minnow enucleated eggs. Theriogenology. 2007b;68: 1282–1291.

57. Arias AM, Brickman JM. Gene expression heterogeneities in embryonic stem cell populations: origin and function. Curr Opin Cell Biol. 2011;23: 650–656. pmid:21982544

58. Latos PA, Sienerth AR, Murray A, Senner CE, Muto M, Ikawa M et al. Elf5-centered transcription factor hub controls trophoblast stem cell self-renewal and differentiation through stoichiometry-sensitive shifts in target gene networks. Genes Dev. 2015;29: 2435–2448. pmid:26584622

59. Muraoka N, Ieda M. Stoichiometry of transcription factors is critical for cardiac reprogramming. Circ Res. 2015;116: 216–218. pmid:25593268

60. Chopra-Dewasthaly R, Korb M, Brunthaler R, Ertl R. Comprehensive RNA-seq profiling to evaluate the sheep mammary gland transcriptome in response to experimental mycoplamsagalactiae infection. PLoS One. 2017;12: E0170015. pmid:28081235

61. Yathish HM, Kumar S, Dubey PP, Modi RP, Chaudhary R, A SK, et al. Profiling of sperm gene transcripts in crossbred (Bos taurus x Bos indicus) bulls. Anim Reprod Sci. 2017;177: 25–34. pmid:27993430

62. Beyhan Z, Lager AE, Cibelli JB. Interspecies nuclear transfer: implications for stem cell biology. Cell Stem Cell. 2007;1: 502–512. pmid:18371390

63. Narbonne P, Miyamoto K, Gurdon JB. Reprogramming and development in nuclear transfer embryos and in interspecific systems. Curr Opin Genet Dev. 2012;22: 450–458. pmid:23062626