Get 20M+ Full-Text Papers For Less Than $1.50/day. Start a 14-Day Trial for You or Your Team.

Learn More →

Identification and characterization of long non-coding RNAs in porcine granulosa cells exposed to 2,3,7,8-tetrachlorodibenzo-p-dioxin

Identification and characterization of long non-coding RNAs in porcine granulosa cells exposed to... Background: Long non-coding RNAs (lncRNAs) may regulate gene expression in numerous biological processes including cellular response to xenobiotics. The exposure of living organisms to 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD), a persistent environmental contaminant, results in reproductive defects in many species including pigs. The aims of the study were to identify and characterize lncRNAs in porcine granulosa cells as well as to examine the effects of TCDD on the lncRNA expression profile in the cells. Results: One thousand six hundred sixty-six lncRNAs were identified and characterized in porcine granulosa cells. The identified lncRNAs were found to be shorter than mRNAs. In addition, the number of exons was lower in lncRNAs than in mRNAs and their exons were longer. TCDD affected the expression of 22 lncRNAs (differentially expressed lncRNAs [DELs]; log fold change ≥ 1, P-adjusted < 0.05) in the examined cells. Potential functions of DELs were indirectly predicted via searching their target cis- and trans-regulated protein-coding genes. The co- expression analysis revealed that DELs may influence the expression of numerous genes, including those involved in cellular response to xenobiotics, dioxin metabolism, endoplasmic reticulum stress and cell proliferation. Aryl hydrocarbon receptor (AhR) and cytochrome P450 1A1 (CYP1A1) were found among the trans-regulated genes. Conclusions: These findings indicate that the identified lncRNAs may constitute a part of the regulatory mechanism of TCDD action in granulosa cells. To our knowledge, this is the first study describing lncRNAs in porcine granulosa cells as well as TCDD effects on the lncRNA expression profile. These results may trigger new research directions leading to better understanding of molecular processes induced by xenobiotics in the ovary. Keywords: AVG-16 cell line, Granulosa cells, lncRNAs, Pig, RNA-Seq, TCDD Background junk or transcriptional noise, but recently they were re- Recent advances in RNA sequencing (RNA-Seq) technolo- ported to play a crucial role in the regulation of gene gies led to the discovery of tens of thousands non-coding expression [1, 2]. The class of ncRNAs includes short RNAs (ncRNAs) across the entire transcriptomes of many ncRNAs (e.g., microRNA, piwiRNA), housekeeping species. Previously, ncRNAs were considered as evolutionary ncRNAs (e.g., ribosomal RNA, transfer RNA) as well as long non-coding RNAs (lncRNAs), which constitute the largest group of ncRNAs [3, 4]. In broad terms, * Correspondence: reniac@uwm.edu.pl 1 lncRNAs are defined as non-coding transcripts longer Department of Animal Anatomy and Physiology, Faculty of Biology and Biotechnology, University of Warmia and Mazury in Olsztyn, Oczapowskiego than 200 nucleotides (nt). Similar to mRNAs, most 1A, 10-719 Olsztyn, Poland known lncRNAs: 1) are transcribed by RNA polymerase Laboratory of Molecular Diagnostics, Faculty of Biology and Biotechnology, II, 2) contain an intron-exon structure and 3) undergo University of Warmia and Mazury in Olsztyn, Prawochenskiego 5, 10-720 Olsztyn, Poland post-transcriptional modifications, e.g., 5′ capping, 3′ Full list of author information is available at the end of the article © The Author(s). 2018 Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated. Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 2 of 13 polyadenylation or alternative splicing [5]. lncRNAs for follicular development, fertilization, implantation originate from intergenic regions of the genome (lincR- and embryo growth [20]. Duetothe AhRpresencein NAs), introns of protein-coding genes (linRNAs) or the porcine granulosa cells [21, 22]andthefactthat TCDD opposite strand of DNA (antisense lncRNAs) [5, 6]. It affects granulosal production of E and P in pigs [14, 2 4 was demonstrated that lncRNA sequences are not 15, 19], it is necessary to indicate molecular targets of highly conserved and their expression profiles differ TCDD in these cells. The effects of TCDD on gene among species, developmental stages and cell types. expression in porcine granulosa cells were described in Due to their variable nature, lncRNAs are still poorly our recent paper [23], but information concerning characterized [7, 8]. lncRNAs in the porcine ovary is scarce [24]. We assumed It was found that lncRNAs regulate gene expression that the transcriptome of porcine granulosa cells contains via numerous mechanisms. lncRNAs may be implicated lncRNAs and that some lncRNAs are involved in TCDD in e.g., chromatin organization, epigenetic modification, action in the cells. Therefore, the aims of the current transcriptional regulation, RNA processing and modifi- study were to: 1) identify and characterize porcine granu- cation, regulation of miRNA activity as well as protein losa lncRNAs, 2) examine the effects of TCDD on the stability or cellular protein localization [3, 8, 9]. lncRNAs lncRNA expression profile (i.e., to identify differentially seem to be involved in many fundamental biological expressed lncRNAs [DELs]) and 3) predict target genes processes including maintenance and induction of stem for DELs in order to tentatively evaluate their regulatory cell pluripotency, cell differentiation, cell proliferation role in granulosal response to TCDD. and apoptosis [6, 8]. The mutations of lncRNAs and deregulation of lncRNA expression have been associated Methods with a number of human diseases, such as different Culture of AVG-16 cells types of cancer as well as neurobiological, cardiovascular AVG-16 cell line originating from granulosa cells of and autoimmune disorders [6, 8]. The results of recent medium porcine follicles was used in the current study studies also indicate that lncRNAs take an active part in (06062701; The European Collection of Authenticated cellular response to various xenobiotics [1]. Cell Cultures; UK; [25]). Previously, we have showed that 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD), a chlori- AVG-16 cells are morphologically and physiologically nated xenobiotic and a member of the polycyclic similar to primary porcine granulosa cells and are a good aromatic hydrocarbon (PAH) family, is a by-product of model for studying xenoestrogen effects on ovarian func- human activity including herbicide production, waste in- tions [22]. Before the experiment, the cells were thawed cineration and fossil fuel burning [10]. Because of TCDD and cultured in six-well plates with seeding density of lipid solubility and high resistance to degradation, its 1×10 cells/3 mL culture medium (Dulbecco’smodi- half-life is long and ranges from 7 to 10 years in human fied Eagle’s medium [DMEM] with 2 mmol/L L-glutamine, bodies and from 25 to 100 years in the environment. As 10% fetal bovine serum [FBS], 0.1 mmol/L non-essential a result, TCDD persists in soil and water sediments as amino acid [NEAA], 2.5 ng/mL fibroblast growth well as in plant and animal organisms [11, 12]. The main factor-basic human [bFGF] and antibiotic mixture [100 U intracellular mechanism of TCDD action involves the penicillin, 100 μg streptomycin and 0.25 μg amphotericin activation of the aryl hydrocarbon receptor (AhR) path- B/mL]; Sigma Aldrich, St. Louis, MO, USA) [22, 23]. After way. AhR is a ligand-activated transcription factor and a reaching 60–65% confluency, the cells were treated with member of the basic helix-loop-helix/PER-ARNT-SIM TCDD (100 nmol/L;Sigma Aldrich) for3,12or24h (bHLH-PAS) family. Upon activation, the ligand-AhR (n = 2 biological replicates per one time point). To complex translocates to the nucleus, dimerizes with AhR bring out the potential of TCDD to transduce intracel- nuclear translocator (ARNT) and induces the expression lular signaling in the examined cells, the selected dose of of different genes including xenobiotic-metabolizing TCDD moderately exceeded its environmentally relevant enzymes [13]. Exposure to TCDD may result in a variety concentration. Moreover, 100 nmol/L of TCDD was found of harmful short- and long-term effects, such as wasting to be effective in porcine granulosa cells (steroidogenesis: syndrome, cancer and neurological dysfunctions. TCDD [26–28]; gene expression: [23, 29]). At the end of the cul- has also been demonstrated to cause endocrine disrup- ture, medium was removed, cells were washed twice with tion and reproductive defects in many species including a phosphate-buffered saline (PBS; Sigma Aldrich) and des- pigs [14–19]. ignated for total RNA isolation. Granulosa cells play a critical role in maintaining ovarian function and, in consequence, female fertility. Total RNA isolation, evaluation of RNA integrity, cDNA These cells protect and nurture oocytes as well as pro- library construction and sequencing duce steroid hormones (estradiol [E ], progesterone Total RNA was isolated from cells using peqGold TriFast [P ]) responsible for creating an optimal environment (Peqlab Biotechnologie GmbH, Erlangen, Germany) 4 Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 3 of 13 according to manufacturer’s instructions. RNA concentra- (score < 0) [38], FlExible Extraction of LncRNAs (FEELnc) tion and quality were determined spectrophotometrically (coding potential > 0.558) [39], Pfam Scan (v1.3) (E-value (NanoVue Plus, GE Healthcare, Little Chalfont, UK). The <0.001) [40] and PLEK (score < 0) [41] to encode proteins integrity of total RNA was assessed by Agilent 2100 Bioa- and 5) transcripts blasted (Blast2GO) to small RNAs nalyzer using RNA 6000 Nano LabChip Kit (Agilent (rRNA, tRNA, snRNA, snoRNA, miRNA, etc.) annotated Technologies, Santa Clara, CA, USA). Only samples with in Rfam database. The remaining transcripts were consid- RNA integrity number (RIN) values higher than 8.0 were ered as lncRNAs and were divided into known lncRNAs used for RNA-Seq. (based on the pig Ensembl database lncRNAs annotation) Depleted RNA obtained from 400 ng of total RNA was and novel lncRNAs. In addition, the obtained lncRNAs used to construct cDNA libraries (TruSeq RNA Sample were classified depending on their genomic locations. The Preparation Kit; Illumina, San Diego, CA, USA). Following protein-coding sequences removed in step 1 of the identi- RNA purification and fragmentation, first and second fication pipeline were characterized and their genomic cDNA strands were synthesized. Next steps included 3’ features were compared with those of lncRNAs (Welch’s end adenylation, adapter ligation and library amplification t-test in R package, P <0.05). (PCR). Quantification of the cDNA library templates was performed using KAPA Library Quantification Kit (Kapa Differential expression analysis of lncRNAs Biosystem, Wilmington, MA, USA). Library profiles were The combined Cufflinks and StringTie tools allowed us estimated with the use of DNA High Sensitivity LabChip to analyze the expression profile of both known and Kit (Agilent Technologies) and 2100 Bioanalyzer. The unannotated lncRNAs. The identified lncRNA sequences cDNA library templates were sequenced using HiSeq 2500 were normalized to FPKM (fragments per kilobase of high throughput sequencing instrument (Illumina) with exon per million fragments mapped) values using the 100 paired-end sequencing. The RNA-Seq analyses of Cuffnorm (version 2.2.1) in the Cufflinks package [32]. protein-coding genes were described in a separate manu- To identify differentially expressed lncRNAs (DELs), the script [23]. This article is focused on the identification and expression levels of lncRNAs in TCDD-treated cells characterization of long non-coding RNAs in porcine were compared to their respective expression levels in granulosa cells exposed to TCDD. control cells. The corresponding P-values were deter- mined for each incubation time by means of Cuffdiff [32]. Sequencing data analysis and transcriptome assembly P-adjusted <0.05 and log fold change (log FC) ≥ 1.0 were 2 2 The sequencing data from this study have been submit- set as a threshold for the significantly different expression. ted to the NCBI BioProject database (http://www.ncb i.nlm.nih.gov/bioproject) under accession number: PRJN Cis- and trans-regulated target gene prediction and Gene A429720. The quality of cDNA fragments obtained after Ontology (GO) enrichment analysis sequencing (raw reads) was first evaluated using FastQC Cis-acting lncRNAs regulate the expression of genes that program (http://www.bioinformatics.babraham.ac.uk/pro are positioned in the vicinity of their transcription sites, jects/fastqc/). Next, the raw reads were trimmed (Trim- whereas trans-acting lncRNAs modulate the expression of momatic tool version 0.32) by removing the adapter genes being at independent loci [42]. Protein-coding genes sequences and reads shorter than 90 nt [30]. After trim- located within 10 kb from DELs were screened with the ming to the same length (90 nt), the fragments were use of R package [43] and selected as DEL potential mapped to the porcine reference genome (Sus_scrofa.Ssc cis-regulated targets [44, 45]. To explore the trans-type rofa10.2; Ensembl database) using STAR version 2.4 [31]. interaction, the Pearson’s correlations coefficients The mapped reads were assembled into contigs with Cuf- (r >0.7 or r < − 0.7) between DELs and protein-coding flinks version 2.2.1 [32] and StringTie version 1.0.4 [33, 34]. genes were analyzed. To analyze functions of the poten- tial DEL target genes and their involvement in biological Identification of lncRNAs processes, the GO database was used (the established cri- A customized multi-step pipeline (Fig. 1a) was employed teria: P-adjusted < 0.05). GO enrichment analysis was per- to identify putative novel lncRNAs in porcine granulosa formed using g:Profiler software [46]. cells. The detailed identification steps included removal of: 1) protein-coding transcripts and transcripts that Alternative splicing events analysis overlap, in sense orientation, with at least one base of all To recognize the differential expression of exons and the porcine protein-coding sequences annotated in the occurrence of splice junctions (i.e., differential usage of Ensembl database; 2) transcripts with a single exon [35, exons and junctions) in the identified lncRNAs, QoRTs 36]; 3) transcripts shorter than 200 nt; 4) transcripts [47] and JunctionSeq [48] tools were applied. First, which were predicted by Coding Potential Calculator QoRTs was used to generate the raw exon- and splice (CPC) (score < 0) [37], Coding-Non-Coding Index (CNCI) junction-counts. Then, JunctionSeq was employed to Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 4 of 13 Fig. 1 Screening of candidate long non-coding RNAs (lncRNAs) in the porcine granulosa cell transcriptome and classification of novel lncRNAs. a Schematic diagram of the pipeline used for the identification of lncRNAs in porcine granulosa cells. b Venn diagram presenting the results of the coding potential analysis of the obtained 4,160 long transcripts. Please note that five different tools (CPC, CNCI, PLEK, Pfam, FEELnc) were employed to analyze the coding potential. In consequence, 1,689 potentially non-coding long transcripts were identified and designated as candidate lncRNAs. c Classification of the obtained novel lncRNAs according to their genomic positions. *Customized algorithms of the authors were applied for the lncRNA identification in this step visualize the identified sites of differential usage of exons synthesized by Thermofisher Scientific Company (TCO as well as splice junction in lncRNA loci of TCDD-treated NS_00031035 – Assay ID: AIPAFN7; TCONS_00016901 – cells in comparison to control cells (P-adjusted <0.05). Assay ID: AIN1HH2; TCONS_00034713 – Assay ID: AIFAT9T). Glyceraldehyde 3-phosphate dehydrogenase Real-time PCR (GAPDH; Assay ID: Ss03373286_u1) and β-actin (Assay ID: To validate RNA-Seq results, three lncRNAs, expression of Ss03376563_uH) were used as reference genes. Real-time which was affected by TCDD in two (TCONS_00016901; PCR was performed using TaqMan® Universal PCR Master TCONS_00031035) or three (TCONS_00034713) incuba- Mix and TaqMan Gene Expression Assay (Thermofisher tion times, were selected for real-time PCR (qRT-PCR). Scientific) in Applied Biosystems 7500 Fast Real-Time PCR Complementary DNA was synthesized using RNA isolated System (Thermofisher Scientific). Thermal cycling condi- from granulosa cells (control and TCDD-treated), 0.5 μmol/ tions consisted of an initial denaturation step at 95 °C for L oligo(dT) primer (Roche, Basel, Switzerland), 1 μmol/L 10 min and then 40 cycles of denaturation at 95 °C for 15 s hexanucleotide primers, 10 U RNase Out (Sigma Aldrich) followed by primer annealing at 60 °C for 1 min. The and Omniscript RT Kit (Qiagen, Hilden, Germany). Reverse qRT-PCR for each sample was carried out in duplicate and transcription reaction was performed at 37 °C for 1 h (Veriti non-template control was included in each run. To present Thermal Cycler, Thermofisher Scientific, Waltham, MA, the data as arbitrary units of the relative expression, lncRNA USA). Specific primers for particular lncRNA were expression levels were normalized to the expression of Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 5 of 13 GAPDH and β-actin. This was done by using comparative applied (Fig. 1a). Twelve thousand one hundred cycle threshold (Ct) method and the quantity based active seventy-three transcripts were obtained after filtering schematic estimating (Q-BASE) model [49]. Data were the protein-coding transcripts and transcripts that expressed as mean ± SEM. The difference in lncRNA overlap with annotated protein-coding sequences. The expression level between control and TCDD-treated cells removal of single exon transcripts and transcripts was evaluated using Student’s t-test (Statistica Software shorter than 200 nt yielded 4,160 transcripts and the Inc., Tulsa, OH, USA). Differences with a probability of evaluation of the protein-coding potential (Fig. 1b) P < 0.05 were considered significant. produced 1,689 potentially non-coding long tran- scripts. The stringency of the selection process was Results additionally increased by discarding transcripts which The sequencing of the porcine granulosa cell transcrip- blasted to small RNAs annotated in the Rfam tome provided 528.5 million reads (50–100 nt). After database. Eventually, 1,666 RNA sequences were discarding adaptor sequences and low-quality sequences identified as lncRNAs, and 42 of them have already (Phred score Q < 20), the remaining 476.2 million reads been annotated in databases. According to genomic were mapped to the annotated whole porcine genome localization of lncRNAs, the identified porcine novel (Sus_scrofa.Sscrofa10.2). The percentage of unique lncRNAs were classified as intergenic lncRNAs (1,319 mapped reads ranged from 76.8% to 80.1%. As a result, transcripts) and antisense lncRNAs (305 transcripts). 56,477 transcripts were collected (see Additional file 1). Thegroup ofantisenselncRNAs (Fig. 1c) were located at the opposite strand of annotated protein-coding Identification of lncRNAs in the porcine granulosa cell genes (281 lncRNAs), snoRNAs (12 lncRNAs), transcriptome snRNAs (7 lncRNAs), miRNAs (3 lncRNAs), To identify lncRNAs from 56,477 mapped transcripts, processed transcript (1 lncRNA) and pseudogene the customized multi-step identification pipeline was (1 lncRNA). Fig. 2 The percentage distribution of the identified lncRNAs and mRNAs according to their (a) transcript length, (b)exon length, (c)exon numberand (d) expression level (the latter expressed as fragments per kilobase of exon per million fragments mapped [FPKM] values). The assumed nucleotide (nt) ranges as well as ranges of FPKM values are depicted in the X axis Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 6 of 13 The comparison of lncRNAs and mRNAs of the granulosa Differentially expressed lncRNAs in TCDD-treated porcine cell transcriptome granulosa cells The transcript length, exon length, exon number and Twenty two differentially expressed lncRNAs (P-adjusted expression level were compared between the identified < 0.05 and log FC ≥ 1.0) were identified. The expression of lncRNAs (1,666 transcripts) and mRNAs (42,913 15, 4 and 7 lncRNAs were significantly altered after 3, 12 transcripts). The length of most of the lncRNAs and 24 h of TCDD treatment, respectively. The log FC ranged from 200 to 1,000 nt (Fig. 2a), while the value for DELs ranged from 3.77 (TCONS_00041714) to length of most lncRNA exons ranged from 50 to − 1.01 (TCONS_00048132) (Table 1). Among all DELs: 1) 200 nt (Fig. 2b). Moreover, a majority of lncRNAs one lncRNA (TCONS_00034713) was up-regulated by consisted of two exons (Fig. 2c)and theexpression TCDD in all incubation times, 2) two lncRNAs (TCO level of more than 50% of lncRNAs was lower than NS_00016901 and TCONS_00031035) were up-regulated FPKM value of 2 (Fig. 2d). by the dioxin after two incubation times (3, 24 and 12, The average length of lncRNAs (832.6 ± 33.06 nt) was 24 h of TCDD treatment, respectively), 3) two lncRNAs shorter (P < 0.05) than that of mRNA (3284.4 ± 12.49 nt) (TCONS_00048132; TCONS_00030731) were down-regu- (Fig. 3a) and the average exon length of lncRNA (318.2 lated, and the downregulation was found after 24 h of ± 10.39 nt) was longer (P < 0.05) than that of mRNA TCDD treatment, 4) five lncRNAs (TCONS_00005658; (295.3 ± 1.03 nt) (Fig. 3b). In addition, the mean exon TCONS_00016901; TCONS_00048979; TCONS_000602 number of lncRNAs (2.6 ± 0.02) was lower (P < 0.05) 23; TCONS_00064401) were expressed only in TCDD- than the mean exon number of mRNAs (11.09 ± 0.05) treated cells and 5) three lncRNAs (TCONS_00038918, (Fig. 3c). The mean expression level (0.51 ± 0.003 vs. TCONS_00030731, TCONS_00064964) up-regulated by 0.49 ± 0.001 FPKM) did not differ (P > 0.05) between TCDD were located in the antisense strand of lncRNA and mRNA (Fig. 3d). protein-coding transcripts (MAF bZIP transcription factor Fig. 3 The comparison of genomic features of the identified lncRNAs and mRNAs. The lncRNAs and mRNAs were compared in respect to average (a) transcript length, (b) exon length, (c) exon number and (d) expression level. Data are expressed as mean ± SE. Statistical analysis was performed using Welch’s t-test in R package. Asterisks designate statistical differences (P < 0.05) Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 7 of 13 Table 1 Differentially expressed lncRNAs (log FC ≥ 1) identified in porcine granulosa cells after 3, 12 or 24 h of TCDD treatment No. Identified lncRNA lncRNA XLOC lncRNA locus log FC 3h 12 h 24h 1 TCONS_00005658 XLOC_002986 10:777270–897286 Inf –– 2 TCONS_00006290 XLOC_003288 7:63840034–63853541 3.34 –– 3 TCONS_00006291 XLOC_003288 10:78713507–78719836 3.32 –– 4 TCONS_00016901 XLOC_008756 14:120019330–120022623 Inf – Inf 5 TCONS_00034713 XLOC_017807 3:108221745–108230820 1.76 1.64 2.32 6 TCONS_00034978 XLOC_017941 3:138550166–138568232 1.16 –– 7 TCONS_00038918 XLOC_019982 5:6915996–7133894 1.34 –– 8 TCONS_00040607 XLOC_020835 5:23100502–23101339 2.07 –– 9 TCONS_00041714 XLOC_021448 6:3237107–3238658 3.77 –– 10 TCONS_00048979 XLOC_024965 7:3272029–3276056 Inf –– 11 TCONS_00053619 XLOC_027505 9:32511923–32566577 1.51 –– 12 TCONS_00056190 XLOC_028876 9:146755972–146762016 2.08 –– 13 TCONS_00060223 XLOC_031885 GL895718.2:3829–4611 Inf –– 14 TCONS_00064401 XLOC_034330 X:123166327–123167560 Inf –– 15 TCONS_00064964 XLOC_034502 X:14951320–15305754 3.40 –– 16 TCONS_00007818 XLOC_004159 11:18740477–18860373 – 1.09 – 17 TCONS_00020891 XLOC_010641 15:71337144–71348412 – 0.99 – 18 TCONS_00031035 XLOC_015924 2:160084074–160102520 – 1.34 2.68 19 TCONS_00008517 XLOC_004581 12:6857198–6863409 –– 2.06 20 TCONS_00030731 XLOC_015756 2:143856229–143907415 –– −0.99 21 TCONS_00031038 XLOC_015924 2:160084074–160102520 –– 2.31 22 TCONS_00048132 XLOC_024540 7:63840034–63853541 –– −1.01 log FC – log fold change 2 2 Inf – transcript expression was not detected in control (untreated) cells F[MafF], transforming growth factor beta-induced to stimulus” (GO:0050896) and “negative regulation of bio- [TGFβI], NHS actin remodeling regulator [NHS]). The logical process” (GO:0048519). The negatively co-expressed expression profile of up- and down-regulated lncRNAs trans-target genes involved, among others, AhR, TGFβI, identified in granulosa cells after 3, 12 or 24 h of TCDD endoplasmic reticulum to nucleus signaling 1 (ERN1), LIM treatment is presented in Fig. 4. domain kinase 2 (LIMK 2) and ephrin type-B receptor 4 (EPHB4) (see Additional file 5). The positively co-expressed The cis- and trans-target genes for DELs genes were enriched in five GO terms (four in “biological The potential target genes for DELs acting in a cis-or process” and one in “molecular function”). Some of the trans-regulatory manner were predicted to investigate genes were related to cellular response to xenobiotics, the possible significance of these lncRNAs in the porcine including “cellular response to chemical stimulus” granulosal response to TCDD. In silico analysis pro- (GO:0070887) and “regulation of signal transduction” duced 10 cis-target genes for DELs (see Additional file 2) (GO:0009966). This group of trans-target genes involved and none of them were enriched (P > 0.05) in the GO cytochrome P4501A1 (CYP1A1), MafF and NHS (see classification. Additional file 6). To analyze the trans-type interaction between DELs (Table 1) and their target genes, the co-expression analysis Alternative splicing events of lncRNAs in TCDD-treated was performed. As a result, a total of 916 negative (see porcine granulosa cells Additional file 3) and 1,143 positive (see Additional file 4) We identified 33 events (sites) of differential usage of exons correlations were detected. The negatively co-expressed and splice junctions in lncRNAs loci of porcine granulosa genes were enriched in 28 GO terms (26 under “biological cells after 3 h of TCDD treatment and only four events process” and two under “molecular function”) including after 12 h of the treatment. No differential usage of exons “intracellular signal transduction” (GO:0035556), “response and junctions was found after 24 h of TCDD treatment. Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 8 of 13 Fig. 4 Heatmap illustrating the expression profile of differentially expressed lncRNAs (presented as Z-score values) in porcine granulosa cells treated with TCDD for 3, 12 or 24 h. Red blocks represent up- and green blocks represent down-regulated lncRNAs. The color scale of the heatmap shows the expression level where the brightest green stands for − 2.0 Z-score and the brightest red stands for + 2.0 Z-score. Z-score was calculated using logFPKM and the following equation: [x -mean(x)] / sd(x)where ‘sd’ means standard deviation. C 3–24: untreated porcine granulosa cells (control) cultured for3, 12or24 h.TCDD 3–24: porcine granulosa cells treated with TCDD (100 nmol/L) for 3, 12 or 24 h. *transcript expression was not detected in control (untreated) cells; ○ transcripts with expression level altered in two or three examined incubation times Among the identified events, 28 and 2 were defined as dif- porcine granulosa cells exposed to TCDD. The effects of ferentially expressed exons, whereas 5 and 2 were described TCDD on the porcine granulosa cell transcriptome have as splice junctions after 3 h and 12 h of TCDD treatment, been previously examined, however the studies have been respectively (see Additional file 7). The JunctionSet profile focused solely on the expression of genes (RNA-Seq [23]; plot for the exemplary lncRNA XLOC (XLOC_020835) is real-time PCR [29]). In comparison to humans and presented in Fig. 5. Two lncRNAs - TCONS_00040606 and rodents, information concerning the identification and TCONS_00040607 – maybeformedonthebasisof this expression of lncRNAs in pigs is limited [24, 50]. To the selected XLOC. Moreover, the expression of lncRNA best of our knowledge, lncRNAs have not yet been identi- TCONS_00040607 correlated negatively with the expres- fied in porcine granulosa cells and TCDD effects on the sion of AhR. lncRNA expression profile have also not been examined. A total of 1,666 lncRNAs and 42,913 mRNAs were Validation of RNA-Seq data by real-time PCR identified in porcine granulosa cells in the present study. To validate the RNA-Seq results, three up-regulated The majority of the identified lncRNAs were classified lncRNAs, i.e., TCONS_00016901, TCONS_00031035 and as intergenic lncRNAs. The identified lncRNAs were TCONS_00034713 were chosen for real-time PCR. The found to be shorter in comparison to those of mRNAs. expression of the selected lncRNAs entirely confirmed the In addition, the exon number of lncRNAs was lower and results obtained by RNA-Seq (Fig. 6). the exon length was longer than those of mRNAs. The features of the identified lncRNAs were similar to those Discussion found in other studies [24, 43, 51, 52]. Recent advances in sequencing technologies revealed The analysis of the expression profile of lncRNAs in that the genome is widely transcribed, yielding a large TCDD-treated porcine granulosa cells revealed the presence number of ncRNAs. In the current study, we investi- of 22 DELs. Similarly to previously reported TCDD-induced gated the genome-wide lncRNA expression profile of changes in granulosal gene expression [23], the most of Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 9 of 13 Fig. 5 Presentation of differential expression of exons (E) and the occurrence of splice junction (J) in the exemplary selected lncRNA XLOC (XLOC_020835) identified in porcine granulosa cells treated with TCDD. The upper panel shows the expression level estimates for the mean normalized read counts for each exon or splice junction of XLOC_020835 identified in TCDD-treated (blue) and control (red) cells. The lower panel shows the exonic regions (boxes, labelled E001-E003) and known splice junctions (solid line, labelled J004). Statistically significant differences (P-adjusted < 0.05) of usage exons and junctions are drawn in violet. C3: untreated porcine granulosa cells (control) cultured for 3 h. TCDD3: porcine granulosa cells treated with TCDD (100 nmol/L) for 3 h DELs (15 out of 22) were detected after 3 h of TCDD results of negatively co-expressed genes suggested a treatment. Since TCDD was shown to induce oxidative possible involvement of DELs in a variety of biological stress in the ovary [53, 54], the altered expression of processes such as “intracellular signal transduction”, lncRNAs after 3 h of TCDD treatment may be associ- “response to stimulus” and “negative regulation of bio- ated with the regulation of the early cellular stress re- logical process”. Interestingly, AhR was found among sponse to TCDD. these negatively trans-regulated genes. Although the Potential functions of DELs identified in porcine gran- downregulation of AhR is a typical response to TCDD, ulosa cells were indirectly predicted via searching their recent reports concerning ncRNAs may provide more target cis- and trans-regulated protein-coding genes. details on TCDD mechanism of action. Li et al. [55] Screening of genes which were located within 10 kb demonstrated that miR-203 (microRNA) suppressed the from DELs enabled us to select 10 cis-target genes, but expression of AhR in TCDD-treated human lung and none of them was enriched into GO classification. To hepatic cells. In the current study, we identified identify the trans-target genes presumptively regulated TCONS_00040607 i.e., lncRNA, the expression of which by lncRNAs, we correlated the expression levels of DELs correlated negatively with that of AhR after 3 h of TCDD and protein-coding genes. A total of 916 negative and treatment. It is possible that TCONS_00040607 affects the 1,143 positive correlations were detected. The GO expression of AhR and is involved in its negative Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 10 of 13 Fig. 6 Real-time PCR validation of three selected differentially expressed lncRNAs (TCONS_00016901, TCONS_00031035, TCONS_00034713) which were identified in TCDD-treated porcine granulosa cells (treated vs. untreated cells) by RNA-Seq. Data are expressed as mean ± SEM (n =3–4). Statistical analysis was performed using Student’s t-test. Asterisks designate statistical differences (P < 0.05). AU: arbitrary units; C: control; TCDD: 2,3,7,8-tetrachlorodibenzo-p-dioxin regulation during the cellular response to TCDD. These playing an important role in the metabolism of xenobi- results suggest that the regulation of AhR by lncRNAs otics [56], correlated positively with the expression of two may constitute part of cellular defense mechanism against DELs (TCONS_00034713 and TCONS_00031305). Due dioxins. This hypothesis, however, needs additional experi- to the fact that TCDD treatment usually induces CYP1A1 mental verification. expression, this enzyme is considered to be a molecular Some of the positively co-expressed genes were enriched marker of TCDD action [57, 58]. In our previous study, in GO terms associated with “cellular response to chem- the expression of CYP1A1 increased significantly in a ical stimulus” and “regulation of signal transduction”.We time-dependent manner after 3, 12 and 24 h of porcine found that the expression of CYP1A1, coding a protein granulosa cell incubation with TCDD [23]. We also Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 11 of 13 demonstrated that the TCDD binding to the porcine events of differential usage of exons and splice junctions in CYP1A1 active site resulted in a rapid closure of the en- lncRNAs loci were identified after 3 h and 12 h of cell zyme substrate channels. This phenomenon may partially incubation with TCDD. The majority of the identified explain TCDD’s high resistance to biodegradation [59]. If events were defined as differentially expressed exons the TCDD binding causes a continuous CYP1A1 blockage, and only a few events were described as splice junc- the cellular response to TCDD may induce an extended tions. It is of interest that two alternatively spliced synthesis of CYP1A1. The fact that the expression of two forms for XLOC_020835 were detected after TCDD DELs: TCONS_00034713 and TCONS_00031305 corre- treatment. In the current study, the expression of one lated positively with the expression of CYP1A1 indicates of these forms i.e., TCONS_00040607, was found to their supportive role in the cellular reaction to TCDD. negatively correlate with the expression of AhR.These Five DELs (TCONS_00005658; TCONS_00016901; facts implicate that TCONS_00040607 is involved in TCONS_00048979; TCONS_00060223; TCONS_0006 the regulation of the TCDD-affected AhR expression in 4401) were found to be expressed only in TCDD-treated porcine granulosa cells. cells. The expression of two DELs (TCONS_00048979 and TCONS_00060223) correlated negatively with the Conclusions expression of the same three genes: ERN1, LIMK2 and In the current study, we identified and characterized EPHB4. ERN1 is associated with endoplasmic reticulum lncRNAs of porcine granulosa cells. We also examined stress and ER protein folding [60, 61]. EPHB4, in turn, is the effects of TCDD on the lncRNA expression profile. linked with the proliferation of ovarian carcinoma cells The co-expression analysis revealed that the identified [62] and LIMK2 with actin microfilament disruption in lncRNAs may influence the expression of numerous porcine oocytes [63]. The obtained data imply that genes, including those involved in: 1) cellular response TCONS_00048979 and TCONS_00060223 are mediators to TCDD (e.g., AhR), 2) dioxin metabolism (e.g., of TCDD action in porcine granulosa cells. CYP1A1, MafF), 3) endoplasmic reticulum stress (e.g., In addition, we identified three DELs located in the ERN1) as well as 4) adhesion and proliferation of cells antisense strand of protein coding genes. These anti- (TGFβ1, NHS, LIMK2). Our results suggest that sense strands were found to positively or negatively lncRNAs constitute a part of the regulatory apparatus of regulate the expression of their sense counterparts and, TCDD action in porcine granulosa cells. This study may therefore, they attracted a lot of attention [64, 65]. In the provide the foundation for future research focused on current study, TCONS_00038918, TCONS_00030731 molecular effects exerted by TCDD in ovarian cells. and TCONS_00064964 were found to be located in the respective antisense strands of MafF, TGFβI and NHS. Additional files MafF belongs to the small MAF family of basic-leucin zipper transcription factors, which affect genes encoding Additional file 1: Summary of the mapping of the RNA sequences to the reference genome. (XLSX 13 kb) proteins responsible for xenobiotic metabolism and Additional file 2: Cis-target genes predicted to be regulated by antioxidation [66–68]. TGFβ1 encodes an extracellular differentially expressed lncRNAs in porcine granulosa cells exposed to matrix (ECM) protein reported to interact with various TCDD for 3, 12 and 24 h. (XLSX 50 kb) matrix molecules (collagens, fibronectin and laminin), Additional file 3: Negatively co-expressed trans-target genes predicted contributing to cell adhesion, proliferation and migration to be regulated by differentially expressed lncRNAs identified in porcine granulosa cells exposed to TCDD for 3, 12 and 24 h. (XLSX 78 kb) [69]. ECM proteins were found to be affected by TCDD Additional file 4: Positively co-expressed trans-target genes predicted in marmosets [70]. NHS products, in turn, were described to be regulated by differentially expressed lncRNAs identified in porcine to maintain cell morphology by remodeling the actin cyto- granulosa cells exposed to TCDD for 3, 12 and 24 h. (XLSX 126 kb) skeleton [71]. The in silico data concerning the possible Additional file 5: Functional enrichment analysis of the negatively relationships between lncRNAs and protein-coding genes co-expresed trans-target genes predicted to be regulated by differentially expressed lncRNAs identified in porcine granulosa cells exposed to TCDD in porcine granulosa cells treated with TCDD are sup- for 3, 12 and 24 h. (XLSX 54 kb) ported by the results of our recent study [72]. In this Additional file 6: Functional enrichment analysis of the positively study, the abundance of heat shock proteins as well as co-expresed trans-target genes predicted to be regulated by differentially cytoskeleton and ECM proteins were significantly affected expressed lncRNAs identified in porcine granulosa cells exposed to TCDD for 3, 12 and 24 h. (XLSX 36 kb) by TCDD in porcine granulosa cells. Additional file 7: Differential expression of exons and the occurence of Similarly to protein-coding genes, lncRNAs also undergo splice junction in XLOCs of the lncRNAs identified in porcine granulosa alternative splicing, resulting in the formation of numerous cells exposed to TCDD for 3 and 12 h. (XLSX 16 kb) lncRNA isoforms and, in consequence, extending their regulatory capabilities [73, 74]. It was demonstrated that Abbreviations xenobiotics may affect alternative splicing, modifying the AhR: Aryl hydrocarbon receptor; antisense lncRNAs: LncRNAs originate from process in a specific manner [75]. In the current study, the the opposite strand of DNA; ARNT: AhR nuclear translocator; bFGF: Fibroblast Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 12 of 13 growth factor-basic human; bHLH-PAS: The basic helix-loop-helix/PER-ARNT- 4. Liu D, Mewalal R, Hu R, Tuskan GA, Yang X. New technologies SIM; CNCI: Coding-Non-Coding Index; CPC: Coding Potential Calculator; accelerate the exploration of non-coding RNAs in horticultural plants. CYP1A1: Cytochrome P450 1A1; DELs: Differentially expressed lncRNAs; Hortic Res. 2017;4:17031. DMEM: Dulbecco’s modified Eagle’s medium; E : Estradiol; EPHB4: Ephrin 5. Prabhakar B, Zhong XB, Rasmussen TP. Exploiting long noncoding RNAs as type-B receptor 4; ERN1: Endoplasmic reticulum to nucleus signaling 1; pharmacological targets to modulate epigenetic diseases. Yale J Biol Med. FBS: Fetal bovine serum; FEELnc: FlExible Extraction of LncRNAs; GO: Gene 2017;90:73–86. Ontology database; LIMK 2: LIM domain kinase 2; lincRNAs: LncRNAs 6. Chen X, Yan CC, Zhang X, You ZH. Long non-coding RNAs and complex originate from intergenic regions of the genome; linRNAs: LncRNAs originate diseases: from experimental results to computational models. Brief from introns of protein-coding genes; lncRNAs: Long non-coding RNAs; Bioinform. 2016;18:558–76. MafF: MAF bZIP transcription factor F; ncRNAs: Non-coding RNAs; 7. Perry RB, Ulitsky I. The functions of long noncoding RNAs in development NEAA: MEM non-essential amino acid solution; NHS: NHS actin remodeling and stem cells. Development. 2016;143:3882–94. regulator; nt: Nucleotides; P : Progesterone; PAH: Polycyclic aromatic 4 8. Szcześniak MW, Makałowska I. lncRNA-RNA interactions across the human hydrocarbon; RIN: RNA integrity number; RNA-Seq: RNA sequencing; Transcriptome. PLoS One. 2016;11(3):e0150353. TCDD: 2,3,7,8-tetrachlorodibenzo-p-dioxin; TGFβI: Transforming growth factor 9. Guo Q, Cheng Y, Liang T, He Y, Ren C, Sun L, et al. Comprehensive beta-induced analysis of lncRNA-mRNA co-expression patterns identifies immune- associated lncRNA biomarkers in ovarian cancer malignant progression. Sci Rep. 2015;5:17683. Acknowledgements 10. Basham KJ, Leonard CJ, Kieffer C, Shelton DN, McDowell ME, Bhonde VR, We would like to thank Piotr Ciereszko for English language review. et al. Dioxin exposure blocks lactation through a direct effect on mammary epithelial cells mediated by the aryl hydrocarbon receptor repressor. Toxicol Funding Sci. 2015;14:36–45. This study was supported by National Science Centre (2012/05/B/NZ9/03333) 11. ATSDR. Draft Update Toxicological Profile for Chlorinated Dibenzo-p-dioxins. and The Ministry of Science and Higher Education in Poland (UWM No. Prepared by Research Triangle Institute for U.S. Department of Health and 528.0206.0806). HumanServices, Agency for Toxic Substances Disease Registry (ATSDR), Atlanta, GA. 677 pp +appendices (1998). Availability of data and materials 12. Nicolopoulou-Stamati P, Pitsos MA. The impact of endocrine disruptors on The datasets analyzed during the current study are available in the NCBI the female reproductive system. Hum Reprod. 2001;7:323–30. BioProject database under accession number: PRJNA429720 (https://www.n 13. Pocar P, Fischer B, Klonisch T, Hombach-Klonisch S. Molecular interactions of cbi.nlm.nih.gov/bioproject/?term=PRJNA429720). the aryl hydrocarbon receptor and its biological and toxicological relevance for reproduction. Reproduction. 2005;129:379–89. Authors’ contributions 14. Piekło R, Grochowalski A, Gregoraszczuk EL. 2,3,7,8-tetrachlorodibenzo- REC, AN, SS & LP contributed to the experimental concept; REC, MR, AN, SS, pdioxin alters follicular steroidogenesis in time- and cell-specific manner. AS & KO participated in creating the experimental design; MR, AS, AN & KO Exp Clin Endocr Diab. 2000;108:299–304. performed cell culture experiments and molecular experiments; LP & JPJ 15. Grochowalski A, Chrzaszcz R, Piekło R, Gregoraszczuk EL. Estrogenic and performed bioinformatics analysis of the data and participated in discussion antiestrogenic effect of in vitro treatment of follicular cells with 2,3,7,8- and interpretation of bioinformatics data; MR, AN & REC analyzed the data, tetrachlorodibenzo-p-dioxin. Chemosphere. 2001;43:823–7. discussed the obtained results and wrote manuscript; AS, KO, SS & TM 16. Petroff BK, Roby KF, Gao X, Son D, Williams S, Johnson D, et al. A participated in the interpretation of the data; MR, REC, AN & JPJ worked out review of mechanisms controlling ovulation with implications for the the figures and tables. All authors read and approved the manuscript. anovulatory effects of polychlorinated dibenzo-p-dioxin in rodents. Toxicology. 2001;158:91–107. 17. Moran FM, Vandevoort CA, Overstreet JW, Lasley BL, Conley AJ. Molecular Ethics approval and consent to participate target of endocrine disruption in human luteinizing granulosa cells by Not applicable. 2,3,7,8-tetrachlorodibenzo-p-dioxin: inhibition of estradiol secretion due to decrease 17a-hydroxylase/17,20-lyase cytochrome P450 expression. Consent for publication Endocrinology. 2003;144:467–73. Not applicable. 18. Mandal PK. Dioxin: a review of its environmental effects and its aryl hydrocarbon receptor biology. J Comp Physiol B. 2005;175:221–30. Competing interests 19. Jablonska O,Piasecka J,Petroff BK,Nynca A,Siawrys G, Wąsowska B, The authors declare that they have no competing interests. et al. In vitro effects of 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) on ovarian, pituitary, and pineal function in pigs. Theriogenology. 2011;76: Author details 921–32. Department of Animal Anatomy and Physiology, Faculty of Biology and 20. Albertini DF, Combelles CM, Benecchi E, Carabatsos MJ. Cellular basis for paracrine Biotechnology, University of Warmia and Mazury in Olsztyn, Oczapowskiego regulation of ovarian follicle development. Reproduction. 2001;121:647–53. 1A, 10-719 Olsztyn, Poland. Laboratory of Molecular Diagnostics, Faculty of 21. Jablonska O, Ciereszko RE. The expression of aryl hydrocarbon receptor in Biology and Biotechnology, University of Warmia and Mazury in Olsztyn, porcine ovarian cells. Reprod Domest Anim. 2013;48:710–6. Prawochenskiego 5, 10-720 Olsztyn, Poland. Department of Plant 22. Sadowska A, Nynca A, Korzeniewska M, Piasecka-Srader J, Jablonska M, Physiology, Genetics and Biotechnology, Faculty of Biology and Orlowska K, et al. Characterization of porcine Granulosa cell line AVG-16. Biotechnology, University of Warmia and Mazury in Olsztyn, Oczapowskiego Folia Biol (Praha). 2015;61:184–94. 1A, 10-719 Olsztyn, Poland. 23. Sadowska A, Nynca A, Ruszkowska M, Paukszto L, Myszczynski K, Orlowska K, et al. Transcriptional profiling of porcine granulosa cells exposed to 2,3,7,8- Received: 7 May 2018 Accepted: 23 August 2018 tetrachlorodibenzo-p-dioxin. Chemosphere. 2017;178:368–77. 24. Tang Z, Wu Y, Yang Y, Yang YT, Wang Z, Yuan J, et al. Comprehensive analysis of long non-coding RNAs highlights their spatio-temporal expression patterns and evolutional conservation in Sus scrofa. Sci Rep. References 2017;7:43166. 1. Dempsey JL, Cui JY. Long non-coding RNAs: a novel paradigm for toxicology. Toxicol Sci. 2017;155:3–21. 25. Horisberger MA. A method for prolonged survival of primary cell lines. In 2. Yu L, Tai L, Zhang L, Chu Y, Li Y, Zhou L. Comparative analyses of long non- Vitro Cell Dev Biol Anim. 2006;42:143–8. coding RNA in lean and obese pigs. Oncotarget. 2017;8:41440–50. 26. Gregoraszczuk EL, Wojtowicz AK, Zabielny E, Grochowalski A. Dose-and time 3. Dykes IM, Emanueli C. Transcriptional and post-transcriptional gene dependent effect of 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) on regulation by long non-coding RNA. Genomics Proteomics Bioinformatics. progesterone secretion by porcine luteal cells cultured in vitro. J Physiol 2017;5:177–86. Pharmacol. 2000;51:127–35. Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 13 of 13 27. Gregoraszczuk EL. Dioxin exposure and porcine reproductive hormonal 52. Xia J, Xin L, Zhu W, Li L, Li C, Wang Y, et al. Characterization of long activity. Cad Saude Publ. 2002;18:453–62. noncoding RNA transcriptome in highenergy diet induced non-alcoholic 28. Jablonska O, Piasecka-Srader J, Nynca A, Kołomycka A, Robak A, Wąsowska steatohepatitis minipigs. Sci Rep. 2016;6:30709. B, et al. 2,3,7,8-Tetrachlorodibenzo-p-dioxin alters steroid secretion but does 53. Melekoglu R, Ciftci O, Cetin A, Basak N, Celik E. The beneficial effects of not affect cell viability and the incidence of apoptosis in porcine luteinised Montelukast against 2,3,7,8-tetrachlorodibenzo-p-dioxin toxicity in female granulosa cells. Acta Vet Hung. 2014;62:408–21. reproductive system in rats. Acta Cir Bras. 2016;31:557–63. 29. Piasecka-Srader J, Sadowska A, Nynca A, Orlowska K, Jablonska M, Jablonska 54. Bhattacharya P, Keating AF. Impact of environmental exposures on ovarian O, et al. The combined effects of 2,3,7,8-tetrachlorodibenzo-p-dioxin and function and role of xenobiotic metabolism during ovotoxicity. Toxicol Appl the phytoestrogen genistein on steroid hormone secretion, AhR and ERβ Pharmacol. 2017;261:227–35. expression and the incidence of apoptosis in granulosa cells of medium 55. Li D, Liu C, Yu H, Zeng X, Xing X, Chenet L, et al. AhR is negatively porcine follicles. J Reprod Dev. 2016;62:103–13. regulated by miR-203 in response to TCDD or BaP treatment. Toxicol Res. 2014;3:142–51. 30. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina 56. Androutsopoulos VP, Tsatsakis AM, Spandidos DA. Cytochrome P450 CYP1A1: sequence data. Bioinformatics. 2014;30:2114–20. wider roles in cancer progression and prevention. BMC Cancer. 2009;9:187. 31. Dobin A, Gingeras TR. Mapping RNA-seq reads with STAR. Curr Protoc 57. Whitlock JP Jr. Induction of cytochrome P4501A1. Annu Rev Pharmacol Bioinformatics. 2015;51:11.14.1–19. Toxicol. 1999;39:103–25. 32. Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, et al. Differential 58. Kim S, Dere E, Burgoon LD, Chang CC, Zacharewski TR. Comparative analysis gene and transcript expression analysis of RNA-seq experiments with of AhR-mediated TCDD-elicited gene expression in human liver adult stem TopHat and cufflinks. Nat Protoc. 2012;7:562–78. cells. Toxicol Sci. 2009;112:229–44. 33. Pertea M, Pertea GM, Antonescu CM, Chang TC, Mendell JT, Salzberg SL. 59. Molcan T, Swigonska S, Orlowska K, Myszczynski K, Nynca A, Sadowska A, StringTie enables improved reconstruction of a transcriptome from RNA-seq et al. Structural-functional adaptations of porcine CYP1A1 to metabolize reads. Nat Biotechnol. 2015;33:290–5. polychlorinated dibenzo-p-dioxins. Chemosphere. 2017;168:205–16. 34. Pertea M, Kim D, Pertea GM, Leek JT, Salzberg SL. Transcript-level expression 60. Minchenko OH. Inhibition of ERN1 signalling enzyme affects hypoxic analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat regulation of the expression of E2F8, EPAS1, HOXC6, ATF3, TBX3 AND Protoc. 2016;11:1650–67. FOXF1 genes in U87 glioma cells. Ukr Biochem J. 2015;87:76–87. 35. Zhan S, Dong Y, Zhao W, Guo J, Zhong T, Wang L, et al. Genome-wide 61. Rashid HO, Yadav RK, Kim HR, Chae HJ. ER stress: autophagy induction, identification and characterization of long non-coding RNAs in inhibition and selection. Autophagy. 2015;11:1956–77. developmental skeletal muscle of fetal goat. BMC Genomics. 2016;7:666. 62. Kumar SR, Masood R, Spannuth WA, Singh J, Scehnet J, Kleiber G, et al. The 36. Han DX, Sun XL, Fu Y, Wang CJ, Liu JB, Jiang H, et al. Identification of long receptor tyrosine kinase EphB4 is overexpressed in ovarian cancer, provides non-coding RNAs in the immature and mature rat anterior pituitary. Sci Rep. survival signals and predicts poor outcome. Br J Cancer. 2007;96:1083–91. 2017;7:17780. 63. Jia RX, Duan X, Song SI, Su SC. LIMK1/2 inhibitor LIMKi 3 suppresses porcine 37. Kong L, Zhang Y, Ye ZQ, Liu XQ, Zhao SQ, Wei L, et al. CPC: assess the oocyte maturation. PeerJ. 2016;4:e2553. protein-coding potential of transcripts using sequence features and support 64. Khorkova O, Myers AJ, Hsiao J, Wahlestedt C. Natural antisense transcripts. vector machine. Nucleic Acids Res. 2007;35:W345–9. Hum Mol Genet. 2014;23(R1):R54–63. 38. Sun L, Luo H, Bu D, Zhao G, Yu K, Zhang C, et al. Utilizing sequence intrinsic 65. Villegas VE, Zaphiropoulos PG. Neighboring gene regulation by antisense composition to classify protein-coding and long non-coding transcripts. long non-coding RNAs. Int J Mol Sci. 2015;16:3251–66. Nucleic Acids Res. 2013;41(17):e166. 66. Dhakshinamoorthy S, Jaiswal AK. Small maf (MafG and MafK) proteins 39. Wucher V, Legeai F, Hédan B, Rizk G, Lagoutte L, Leeb T, et al. FEELnc: a tool negatively regulate antioxidant response element-mediated expression and for long non-coding RNA annotation and its application to the dog antioxidant induction of the NAD(P)H:Quinone oxidoreductase1 gene. J Biol transcriptome. Nucleic Acids Res. 2017;45(8):e57. Chem. 2000;275:40134–41. 40. Mistry J, Bateman A, Finn RD. Predicting active site residue annotations in 67. Katsuoka F, Motohashi H, Ishii T, Aburatani H, Engel JD, Yamamoto M. Genetic the Pfam database. BMC Bioinformatics. 2007;8:298. evidence that small maf proteins are essential for the activation of antioxidant 41. Li A, Zhang J, Zhou Z. PLEK: a tool for predicting long non-coding RNAs response element-dependent genes. Mol Cell Biol. 2005;25:8044–51. and messenger RNAs based on an improved k-mer scheme. BMC 68. Massrieh W, Derjuga A, Doualla-Bell F, Ku CY, Sanborn BM, Blank V. Bioinformatics. 2014;15:311. Regulation of the MAFF transcription factor by proinflammatory cytokines in 42. Kung JT, Colognori D, Lee JT. Long noncoding RNAs: past, present, and myometrial cells. Biol Reprod. 2006;74:699–705. future. Genetics. 2013;193:651–69. 69. Kim JE. Molecular properties of wild-type and mutant betaIG-H3 proteins. 43. Huber W, Carey VJ, Gentleman R, Anders S, Carlson M, Carvalho BS, et al. Invest Ophthalmol Vis Sci. 2002;43:656–61. Orchestrating high-throughput genomic analysis with bioconductor. Nat 70. Riecke K, Grimm D, Shakibaei M, Kossmehl P, Schulze-Tanzil G, Paul M, et al. Methods. 2015;12:115–21. Low doses of 2,3,7,8-tetrachlorodibenzo-p-dioxin increase transforming 44. Ren H, Wang G, Chen L, Jiang J, Liu L, Li N, et al. Genome-wide analysis of growth factor beta and cause myocardial fibrosis in marmosets long non-coding RNAs at early stage of skin pigmentation in goats (Capra (Callithrixjacchus). Arch Toxicol. 2002;76:360–6. hircus). BMC Genomics. 2016;17:67. 71. Zhang DD, Du JZ, Topolewski J, Wang XM. Review recent progress in 45. Zhang T, Zhang X, Han K, Zhang G, Wang J, Xie K, et al. Genome-wide identification and characterization of loci associated with sex-linked analysis of lncRNA and mRNA expression during differentiation of congenital cataract. Genet Mol Res. 2016;29, 15(3) abdominal Preadipocytes in the chicken. G3 (Bethesda). 2017;7:953–66. 72. Orlowska K, Swigonska S, Sadowska A, Ruszkowska M, Nynca A, Molcan T, 46. Reimand J, Kull M, Peterson H, Hansen J, Vilo J. G:profiler—a web-based et al. The effects of 2,3,7,8-tetrachlorodibenzo-p-dioxin on the proteome of toolset for functional profiling of gene lists from large-scale experiments. porcine granulosa cells. Chemosphere. 2018;212:170–81. Nucleic Acids Res. 2007;35:W193–200. 73. Bozgeyik E, Igci YZ, Sami Jacksi MF, Arman K, Gurses SA, Bozgeyik I, et al. A 47. Hartley SW, Mullikin JC. QoRTs: a comprehensive toolset for quality novel variable exonic region and differential expression of LINC00663 non- control and data processing of RNA-Seq experiments. BMC Bioinformatics. coding RNA in various cancer cell lines and normal human tissue samples. 2015;16:224. Tumour Biol. 2016;37:8791–8. 48. Hartley SW, Mullikin JC. Detection and visualization of differential splicing in 74. Sen R, Doose G, Stadler PT. Rare splice variants in long non-coding RNAs. RNA-Seq data with JunctionSeq. Nucleic Acids Res. 2016;44:e127. Non-coding RNA. 2017;3(3):23. 49. Hellemans J,Mortier G, De PaepeA,Speleman F,VandesompeleJ. qBase relative 75. Zaharieva E, Chipman JK, Soller M. Alternative splicing interference by quantification framework and software for management and automated analysis xenobiotics. Toxicology. 2012;296:1–12. of real-time quantitative PCR data. Genome Biol 2007;8(2):R19. 50. Zhao W, Mu Y, Ma L, Wang C, Tang Z, Yang S, et al. Systematic identification and characterization of long intergenic non-coding RNAs in fetal porcine skeletal muscle development. Sci Rep. 2015;5:8957. 51. Ran M, Chen B, Li Z, Wu M, Liu X, He C, et al. Systematic identification of long noncoding RNAs in immature and mature porcine testes. Biol Reprod. 2016;94:77. http://www.deepdyve.com/assets/images/DeepDyve-Logo-lg.png Journal of Animal Science and Biotechnology Springer Journals

Identification and characterization of long non-coding RNAs in porcine granulosa cells exposed to 2,3,7,8-tetrachlorodibenzo-p-dioxin

Download PDF
13 pages

Loading next page...
 
/lp/springer-journals/identification-and-characterization-of-long-non-coding-rnas-in-porcine-jbW5yb2lvH

References (47)

Publisher
Springer Journals
Copyright
Copyright © 2018 by The Author(s).
Subject
Life Sciences; Agriculture; Biotechnology; Food Science; Animal Genetics and Genomics; Animal Physiology
eISSN
2049-1891
DOI
10.1186/s40104-018-0288-3
Publisher site
See Article on Publisher Site

Abstract

Background: Long non-coding RNAs (lncRNAs) may regulate gene expression in numerous biological processes including cellular response to xenobiotics. The exposure of living organisms to 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD), a persistent environmental contaminant, results in reproductive defects in many species including pigs. The aims of the study were to identify and characterize lncRNAs in porcine granulosa cells as well as to examine the effects of TCDD on the lncRNA expression profile in the cells. Results: One thousand six hundred sixty-six lncRNAs were identified and characterized in porcine granulosa cells. The identified lncRNAs were found to be shorter than mRNAs. In addition, the number of exons was lower in lncRNAs than in mRNAs and their exons were longer. TCDD affected the expression of 22 lncRNAs (differentially expressed lncRNAs [DELs]; log fold change ≥ 1, P-adjusted < 0.05) in the examined cells. Potential functions of DELs were indirectly predicted via searching their target cis- and trans-regulated protein-coding genes. The co- expression analysis revealed that DELs may influence the expression of numerous genes, including those involved in cellular response to xenobiotics, dioxin metabolism, endoplasmic reticulum stress and cell proliferation. Aryl hydrocarbon receptor (AhR) and cytochrome P450 1A1 (CYP1A1) were found among the trans-regulated genes. Conclusions: These findings indicate that the identified lncRNAs may constitute a part of the regulatory mechanism of TCDD action in granulosa cells. To our knowledge, this is the first study describing lncRNAs in porcine granulosa cells as well as TCDD effects on the lncRNA expression profile. These results may trigger new research directions leading to better understanding of molecular processes induced by xenobiotics in the ovary. Keywords: AVG-16 cell line, Granulosa cells, lncRNAs, Pig, RNA-Seq, TCDD Background junk or transcriptional noise, but recently they were re- Recent advances in RNA sequencing (RNA-Seq) technolo- ported to play a crucial role in the regulation of gene gies led to the discovery of tens of thousands non-coding expression [1, 2]. The class of ncRNAs includes short RNAs (ncRNAs) across the entire transcriptomes of many ncRNAs (e.g., microRNA, piwiRNA), housekeeping species. Previously, ncRNAs were considered as evolutionary ncRNAs (e.g., ribosomal RNA, transfer RNA) as well as long non-coding RNAs (lncRNAs), which constitute the largest group of ncRNAs [3, 4]. In broad terms, * Correspondence: reniac@uwm.edu.pl 1 lncRNAs are defined as non-coding transcripts longer Department of Animal Anatomy and Physiology, Faculty of Biology and Biotechnology, University of Warmia and Mazury in Olsztyn, Oczapowskiego than 200 nucleotides (nt). Similar to mRNAs, most 1A, 10-719 Olsztyn, Poland known lncRNAs: 1) are transcribed by RNA polymerase Laboratory of Molecular Diagnostics, Faculty of Biology and Biotechnology, II, 2) contain an intron-exon structure and 3) undergo University of Warmia and Mazury in Olsztyn, Prawochenskiego 5, 10-720 Olsztyn, Poland post-transcriptional modifications, e.g., 5′ capping, 3′ Full list of author information is available at the end of the article © The Author(s). 2018 Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated. Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 2 of 13 polyadenylation or alternative splicing [5]. lncRNAs for follicular development, fertilization, implantation originate from intergenic regions of the genome (lincR- and embryo growth [20]. Duetothe AhRpresencein NAs), introns of protein-coding genes (linRNAs) or the porcine granulosa cells [21, 22]andthefactthat TCDD opposite strand of DNA (antisense lncRNAs) [5, 6]. It affects granulosal production of E and P in pigs [14, 2 4 was demonstrated that lncRNA sequences are not 15, 19], it is necessary to indicate molecular targets of highly conserved and their expression profiles differ TCDD in these cells. The effects of TCDD on gene among species, developmental stages and cell types. expression in porcine granulosa cells were described in Due to their variable nature, lncRNAs are still poorly our recent paper [23], but information concerning characterized [7, 8]. lncRNAs in the porcine ovary is scarce [24]. We assumed It was found that lncRNAs regulate gene expression that the transcriptome of porcine granulosa cells contains via numerous mechanisms. lncRNAs may be implicated lncRNAs and that some lncRNAs are involved in TCDD in e.g., chromatin organization, epigenetic modification, action in the cells. Therefore, the aims of the current transcriptional regulation, RNA processing and modifi- study were to: 1) identify and characterize porcine granu- cation, regulation of miRNA activity as well as protein losa lncRNAs, 2) examine the effects of TCDD on the stability or cellular protein localization [3, 8, 9]. lncRNAs lncRNA expression profile (i.e., to identify differentially seem to be involved in many fundamental biological expressed lncRNAs [DELs]) and 3) predict target genes processes including maintenance and induction of stem for DELs in order to tentatively evaluate their regulatory cell pluripotency, cell differentiation, cell proliferation role in granulosal response to TCDD. and apoptosis [6, 8]. The mutations of lncRNAs and deregulation of lncRNA expression have been associated Methods with a number of human diseases, such as different Culture of AVG-16 cells types of cancer as well as neurobiological, cardiovascular AVG-16 cell line originating from granulosa cells of and autoimmune disorders [6, 8]. The results of recent medium porcine follicles was used in the current study studies also indicate that lncRNAs take an active part in (06062701; The European Collection of Authenticated cellular response to various xenobiotics [1]. Cell Cultures; UK; [25]). Previously, we have showed that 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD), a chlori- AVG-16 cells are morphologically and physiologically nated xenobiotic and a member of the polycyclic similar to primary porcine granulosa cells and are a good aromatic hydrocarbon (PAH) family, is a by-product of model for studying xenoestrogen effects on ovarian func- human activity including herbicide production, waste in- tions [22]. Before the experiment, the cells were thawed cineration and fossil fuel burning [10]. Because of TCDD and cultured in six-well plates with seeding density of lipid solubility and high resistance to degradation, its 1×10 cells/3 mL culture medium (Dulbecco’smodi- half-life is long and ranges from 7 to 10 years in human fied Eagle’s medium [DMEM] with 2 mmol/L L-glutamine, bodies and from 25 to 100 years in the environment. As 10% fetal bovine serum [FBS], 0.1 mmol/L non-essential a result, TCDD persists in soil and water sediments as amino acid [NEAA], 2.5 ng/mL fibroblast growth well as in plant and animal organisms [11, 12]. The main factor-basic human [bFGF] and antibiotic mixture [100 U intracellular mechanism of TCDD action involves the penicillin, 100 μg streptomycin and 0.25 μg amphotericin activation of the aryl hydrocarbon receptor (AhR) path- B/mL]; Sigma Aldrich, St. Louis, MO, USA) [22, 23]. After way. AhR is a ligand-activated transcription factor and a reaching 60–65% confluency, the cells were treated with member of the basic helix-loop-helix/PER-ARNT-SIM TCDD (100 nmol/L;Sigma Aldrich) for3,12or24h (bHLH-PAS) family. Upon activation, the ligand-AhR (n = 2 biological replicates per one time point). To complex translocates to the nucleus, dimerizes with AhR bring out the potential of TCDD to transduce intracel- nuclear translocator (ARNT) and induces the expression lular signaling in the examined cells, the selected dose of of different genes including xenobiotic-metabolizing TCDD moderately exceeded its environmentally relevant enzymes [13]. Exposure to TCDD may result in a variety concentration. Moreover, 100 nmol/L of TCDD was found of harmful short- and long-term effects, such as wasting to be effective in porcine granulosa cells (steroidogenesis: syndrome, cancer and neurological dysfunctions. TCDD [26–28]; gene expression: [23, 29]). At the end of the cul- has also been demonstrated to cause endocrine disrup- ture, medium was removed, cells were washed twice with tion and reproductive defects in many species including a phosphate-buffered saline (PBS; Sigma Aldrich) and des- pigs [14–19]. ignated for total RNA isolation. Granulosa cells play a critical role in maintaining ovarian function and, in consequence, female fertility. Total RNA isolation, evaluation of RNA integrity, cDNA These cells protect and nurture oocytes as well as pro- library construction and sequencing duce steroid hormones (estradiol [E ], progesterone Total RNA was isolated from cells using peqGold TriFast [P ]) responsible for creating an optimal environment (Peqlab Biotechnologie GmbH, Erlangen, Germany) 4 Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 3 of 13 according to manufacturer’s instructions. RNA concentra- (score < 0) [38], FlExible Extraction of LncRNAs (FEELnc) tion and quality were determined spectrophotometrically (coding potential > 0.558) [39], Pfam Scan (v1.3) (E-value (NanoVue Plus, GE Healthcare, Little Chalfont, UK). The <0.001) [40] and PLEK (score < 0) [41] to encode proteins integrity of total RNA was assessed by Agilent 2100 Bioa- and 5) transcripts blasted (Blast2GO) to small RNAs nalyzer using RNA 6000 Nano LabChip Kit (Agilent (rRNA, tRNA, snRNA, snoRNA, miRNA, etc.) annotated Technologies, Santa Clara, CA, USA). Only samples with in Rfam database. The remaining transcripts were consid- RNA integrity number (RIN) values higher than 8.0 were ered as lncRNAs and were divided into known lncRNAs used for RNA-Seq. (based on the pig Ensembl database lncRNAs annotation) Depleted RNA obtained from 400 ng of total RNA was and novel lncRNAs. In addition, the obtained lncRNAs used to construct cDNA libraries (TruSeq RNA Sample were classified depending on their genomic locations. The Preparation Kit; Illumina, San Diego, CA, USA). Following protein-coding sequences removed in step 1 of the identi- RNA purification and fragmentation, first and second fication pipeline were characterized and their genomic cDNA strands were synthesized. Next steps included 3’ features were compared with those of lncRNAs (Welch’s end adenylation, adapter ligation and library amplification t-test in R package, P <0.05). (PCR). Quantification of the cDNA library templates was performed using KAPA Library Quantification Kit (Kapa Differential expression analysis of lncRNAs Biosystem, Wilmington, MA, USA). Library profiles were The combined Cufflinks and StringTie tools allowed us estimated with the use of DNA High Sensitivity LabChip to analyze the expression profile of both known and Kit (Agilent Technologies) and 2100 Bioanalyzer. The unannotated lncRNAs. The identified lncRNA sequences cDNA library templates were sequenced using HiSeq 2500 were normalized to FPKM (fragments per kilobase of high throughput sequencing instrument (Illumina) with exon per million fragments mapped) values using the 100 paired-end sequencing. The RNA-Seq analyses of Cuffnorm (version 2.2.1) in the Cufflinks package [32]. protein-coding genes were described in a separate manu- To identify differentially expressed lncRNAs (DELs), the script [23]. This article is focused on the identification and expression levels of lncRNAs in TCDD-treated cells characterization of long non-coding RNAs in porcine were compared to their respective expression levels in granulosa cells exposed to TCDD. control cells. The corresponding P-values were deter- mined for each incubation time by means of Cuffdiff [32]. Sequencing data analysis and transcriptome assembly P-adjusted <0.05 and log fold change (log FC) ≥ 1.0 were 2 2 The sequencing data from this study have been submit- set as a threshold for the significantly different expression. ted to the NCBI BioProject database (http://www.ncb i.nlm.nih.gov/bioproject) under accession number: PRJN Cis- and trans-regulated target gene prediction and Gene A429720. The quality of cDNA fragments obtained after Ontology (GO) enrichment analysis sequencing (raw reads) was first evaluated using FastQC Cis-acting lncRNAs regulate the expression of genes that program (http://www.bioinformatics.babraham.ac.uk/pro are positioned in the vicinity of their transcription sites, jects/fastqc/). Next, the raw reads were trimmed (Trim- whereas trans-acting lncRNAs modulate the expression of momatic tool version 0.32) by removing the adapter genes being at independent loci [42]. Protein-coding genes sequences and reads shorter than 90 nt [30]. After trim- located within 10 kb from DELs were screened with the ming to the same length (90 nt), the fragments were use of R package [43] and selected as DEL potential mapped to the porcine reference genome (Sus_scrofa.Ssc cis-regulated targets [44, 45]. To explore the trans-type rofa10.2; Ensembl database) using STAR version 2.4 [31]. interaction, the Pearson’s correlations coefficients The mapped reads were assembled into contigs with Cuf- (r >0.7 or r < − 0.7) between DELs and protein-coding flinks version 2.2.1 [32] and StringTie version 1.0.4 [33, 34]. genes were analyzed. To analyze functions of the poten- tial DEL target genes and their involvement in biological Identification of lncRNAs processes, the GO database was used (the established cri- A customized multi-step pipeline (Fig. 1a) was employed teria: P-adjusted < 0.05). GO enrichment analysis was per- to identify putative novel lncRNAs in porcine granulosa formed using g:Profiler software [46]. cells. The detailed identification steps included removal of: 1) protein-coding transcripts and transcripts that Alternative splicing events analysis overlap, in sense orientation, with at least one base of all To recognize the differential expression of exons and the porcine protein-coding sequences annotated in the occurrence of splice junctions (i.e., differential usage of Ensembl database; 2) transcripts with a single exon [35, exons and junctions) in the identified lncRNAs, QoRTs 36]; 3) transcripts shorter than 200 nt; 4) transcripts [47] and JunctionSeq [48] tools were applied. First, which were predicted by Coding Potential Calculator QoRTs was used to generate the raw exon- and splice (CPC) (score < 0) [37], Coding-Non-Coding Index (CNCI) junction-counts. Then, JunctionSeq was employed to Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 4 of 13 Fig. 1 Screening of candidate long non-coding RNAs (lncRNAs) in the porcine granulosa cell transcriptome and classification of novel lncRNAs. a Schematic diagram of the pipeline used for the identification of lncRNAs in porcine granulosa cells. b Venn diagram presenting the results of the coding potential analysis of the obtained 4,160 long transcripts. Please note that five different tools (CPC, CNCI, PLEK, Pfam, FEELnc) were employed to analyze the coding potential. In consequence, 1,689 potentially non-coding long transcripts were identified and designated as candidate lncRNAs. c Classification of the obtained novel lncRNAs according to their genomic positions. *Customized algorithms of the authors were applied for the lncRNA identification in this step visualize the identified sites of differential usage of exons synthesized by Thermofisher Scientific Company (TCO as well as splice junction in lncRNA loci of TCDD-treated NS_00031035 – Assay ID: AIPAFN7; TCONS_00016901 – cells in comparison to control cells (P-adjusted <0.05). Assay ID: AIN1HH2; TCONS_00034713 – Assay ID: AIFAT9T). Glyceraldehyde 3-phosphate dehydrogenase Real-time PCR (GAPDH; Assay ID: Ss03373286_u1) and β-actin (Assay ID: To validate RNA-Seq results, three lncRNAs, expression of Ss03376563_uH) were used as reference genes. Real-time which was affected by TCDD in two (TCONS_00016901; PCR was performed using TaqMan® Universal PCR Master TCONS_00031035) or three (TCONS_00034713) incuba- Mix and TaqMan Gene Expression Assay (Thermofisher tion times, were selected for real-time PCR (qRT-PCR). Scientific) in Applied Biosystems 7500 Fast Real-Time PCR Complementary DNA was synthesized using RNA isolated System (Thermofisher Scientific). Thermal cycling condi- from granulosa cells (control and TCDD-treated), 0.5 μmol/ tions consisted of an initial denaturation step at 95 °C for L oligo(dT) primer (Roche, Basel, Switzerland), 1 μmol/L 10 min and then 40 cycles of denaturation at 95 °C for 15 s hexanucleotide primers, 10 U RNase Out (Sigma Aldrich) followed by primer annealing at 60 °C for 1 min. The and Omniscript RT Kit (Qiagen, Hilden, Germany). Reverse qRT-PCR for each sample was carried out in duplicate and transcription reaction was performed at 37 °C for 1 h (Veriti non-template control was included in each run. To present Thermal Cycler, Thermofisher Scientific, Waltham, MA, the data as arbitrary units of the relative expression, lncRNA USA). Specific primers for particular lncRNA were expression levels were normalized to the expression of Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 5 of 13 GAPDH and β-actin. This was done by using comparative applied (Fig. 1a). Twelve thousand one hundred cycle threshold (Ct) method and the quantity based active seventy-three transcripts were obtained after filtering schematic estimating (Q-BASE) model [49]. Data were the protein-coding transcripts and transcripts that expressed as mean ± SEM. The difference in lncRNA overlap with annotated protein-coding sequences. The expression level between control and TCDD-treated cells removal of single exon transcripts and transcripts was evaluated using Student’s t-test (Statistica Software shorter than 200 nt yielded 4,160 transcripts and the Inc., Tulsa, OH, USA). Differences with a probability of evaluation of the protein-coding potential (Fig. 1b) P < 0.05 were considered significant. produced 1,689 potentially non-coding long tran- scripts. The stringency of the selection process was Results additionally increased by discarding transcripts which The sequencing of the porcine granulosa cell transcrip- blasted to small RNAs annotated in the Rfam tome provided 528.5 million reads (50–100 nt). After database. Eventually, 1,666 RNA sequences were discarding adaptor sequences and low-quality sequences identified as lncRNAs, and 42 of them have already (Phred score Q < 20), the remaining 476.2 million reads been annotated in databases. According to genomic were mapped to the annotated whole porcine genome localization of lncRNAs, the identified porcine novel (Sus_scrofa.Sscrofa10.2). The percentage of unique lncRNAs were classified as intergenic lncRNAs (1,319 mapped reads ranged from 76.8% to 80.1%. As a result, transcripts) and antisense lncRNAs (305 transcripts). 56,477 transcripts were collected (see Additional file 1). Thegroup ofantisenselncRNAs (Fig. 1c) were located at the opposite strand of annotated protein-coding Identification of lncRNAs in the porcine granulosa cell genes (281 lncRNAs), snoRNAs (12 lncRNAs), transcriptome snRNAs (7 lncRNAs), miRNAs (3 lncRNAs), To identify lncRNAs from 56,477 mapped transcripts, processed transcript (1 lncRNA) and pseudogene the customized multi-step identification pipeline was (1 lncRNA). Fig. 2 The percentage distribution of the identified lncRNAs and mRNAs according to their (a) transcript length, (b)exon length, (c)exon numberand (d) expression level (the latter expressed as fragments per kilobase of exon per million fragments mapped [FPKM] values). The assumed nucleotide (nt) ranges as well as ranges of FPKM values are depicted in the X axis Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 6 of 13 The comparison of lncRNAs and mRNAs of the granulosa Differentially expressed lncRNAs in TCDD-treated porcine cell transcriptome granulosa cells The transcript length, exon length, exon number and Twenty two differentially expressed lncRNAs (P-adjusted expression level were compared between the identified < 0.05 and log FC ≥ 1.0) were identified. The expression of lncRNAs (1,666 transcripts) and mRNAs (42,913 15, 4 and 7 lncRNAs were significantly altered after 3, 12 transcripts). The length of most of the lncRNAs and 24 h of TCDD treatment, respectively. The log FC ranged from 200 to 1,000 nt (Fig. 2a), while the value for DELs ranged from 3.77 (TCONS_00041714) to length of most lncRNA exons ranged from 50 to − 1.01 (TCONS_00048132) (Table 1). Among all DELs: 1) 200 nt (Fig. 2b). Moreover, a majority of lncRNAs one lncRNA (TCONS_00034713) was up-regulated by consisted of two exons (Fig. 2c)and theexpression TCDD in all incubation times, 2) two lncRNAs (TCO level of more than 50% of lncRNAs was lower than NS_00016901 and TCONS_00031035) were up-regulated FPKM value of 2 (Fig. 2d). by the dioxin after two incubation times (3, 24 and 12, The average length of lncRNAs (832.6 ± 33.06 nt) was 24 h of TCDD treatment, respectively), 3) two lncRNAs shorter (P < 0.05) than that of mRNA (3284.4 ± 12.49 nt) (TCONS_00048132; TCONS_00030731) were down-regu- (Fig. 3a) and the average exon length of lncRNA (318.2 lated, and the downregulation was found after 24 h of ± 10.39 nt) was longer (P < 0.05) than that of mRNA TCDD treatment, 4) five lncRNAs (TCONS_00005658; (295.3 ± 1.03 nt) (Fig. 3b). In addition, the mean exon TCONS_00016901; TCONS_00048979; TCONS_000602 number of lncRNAs (2.6 ± 0.02) was lower (P < 0.05) 23; TCONS_00064401) were expressed only in TCDD- than the mean exon number of mRNAs (11.09 ± 0.05) treated cells and 5) three lncRNAs (TCONS_00038918, (Fig. 3c). The mean expression level (0.51 ± 0.003 vs. TCONS_00030731, TCONS_00064964) up-regulated by 0.49 ± 0.001 FPKM) did not differ (P > 0.05) between TCDD were located in the antisense strand of lncRNA and mRNA (Fig. 3d). protein-coding transcripts (MAF bZIP transcription factor Fig. 3 The comparison of genomic features of the identified lncRNAs and mRNAs. The lncRNAs and mRNAs were compared in respect to average (a) transcript length, (b) exon length, (c) exon number and (d) expression level. Data are expressed as mean ± SE. Statistical analysis was performed using Welch’s t-test in R package. Asterisks designate statistical differences (P < 0.05) Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 7 of 13 Table 1 Differentially expressed lncRNAs (log FC ≥ 1) identified in porcine granulosa cells after 3, 12 or 24 h of TCDD treatment No. Identified lncRNA lncRNA XLOC lncRNA locus log FC 3h 12 h 24h 1 TCONS_00005658 XLOC_002986 10:777270–897286 Inf –– 2 TCONS_00006290 XLOC_003288 7:63840034–63853541 3.34 –– 3 TCONS_00006291 XLOC_003288 10:78713507–78719836 3.32 –– 4 TCONS_00016901 XLOC_008756 14:120019330–120022623 Inf – Inf 5 TCONS_00034713 XLOC_017807 3:108221745–108230820 1.76 1.64 2.32 6 TCONS_00034978 XLOC_017941 3:138550166–138568232 1.16 –– 7 TCONS_00038918 XLOC_019982 5:6915996–7133894 1.34 –– 8 TCONS_00040607 XLOC_020835 5:23100502–23101339 2.07 –– 9 TCONS_00041714 XLOC_021448 6:3237107–3238658 3.77 –– 10 TCONS_00048979 XLOC_024965 7:3272029–3276056 Inf –– 11 TCONS_00053619 XLOC_027505 9:32511923–32566577 1.51 –– 12 TCONS_00056190 XLOC_028876 9:146755972–146762016 2.08 –– 13 TCONS_00060223 XLOC_031885 GL895718.2:3829–4611 Inf –– 14 TCONS_00064401 XLOC_034330 X:123166327–123167560 Inf –– 15 TCONS_00064964 XLOC_034502 X:14951320–15305754 3.40 –– 16 TCONS_00007818 XLOC_004159 11:18740477–18860373 – 1.09 – 17 TCONS_00020891 XLOC_010641 15:71337144–71348412 – 0.99 – 18 TCONS_00031035 XLOC_015924 2:160084074–160102520 – 1.34 2.68 19 TCONS_00008517 XLOC_004581 12:6857198–6863409 –– 2.06 20 TCONS_00030731 XLOC_015756 2:143856229–143907415 –– −0.99 21 TCONS_00031038 XLOC_015924 2:160084074–160102520 –– 2.31 22 TCONS_00048132 XLOC_024540 7:63840034–63853541 –– −1.01 log FC – log fold change 2 2 Inf – transcript expression was not detected in control (untreated) cells F[MafF], transforming growth factor beta-induced to stimulus” (GO:0050896) and “negative regulation of bio- [TGFβI], NHS actin remodeling regulator [NHS]). The logical process” (GO:0048519). The negatively co-expressed expression profile of up- and down-regulated lncRNAs trans-target genes involved, among others, AhR, TGFβI, identified in granulosa cells after 3, 12 or 24 h of TCDD endoplasmic reticulum to nucleus signaling 1 (ERN1), LIM treatment is presented in Fig. 4. domain kinase 2 (LIMK 2) and ephrin type-B receptor 4 (EPHB4) (see Additional file 5). The positively co-expressed The cis- and trans-target genes for DELs genes were enriched in five GO terms (four in “biological The potential target genes for DELs acting in a cis-or process” and one in “molecular function”). Some of the trans-regulatory manner were predicted to investigate genes were related to cellular response to xenobiotics, the possible significance of these lncRNAs in the porcine including “cellular response to chemical stimulus” granulosal response to TCDD. In silico analysis pro- (GO:0070887) and “regulation of signal transduction” duced 10 cis-target genes for DELs (see Additional file 2) (GO:0009966). This group of trans-target genes involved and none of them were enriched (P > 0.05) in the GO cytochrome P4501A1 (CYP1A1), MafF and NHS (see classification. Additional file 6). To analyze the trans-type interaction between DELs (Table 1) and their target genes, the co-expression analysis Alternative splicing events of lncRNAs in TCDD-treated was performed. As a result, a total of 916 negative (see porcine granulosa cells Additional file 3) and 1,143 positive (see Additional file 4) We identified 33 events (sites) of differential usage of exons correlations were detected. The negatively co-expressed and splice junctions in lncRNAs loci of porcine granulosa genes were enriched in 28 GO terms (26 under “biological cells after 3 h of TCDD treatment and only four events process” and two under “molecular function”) including after 12 h of the treatment. No differential usage of exons “intracellular signal transduction” (GO:0035556), “response and junctions was found after 24 h of TCDD treatment. Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 8 of 13 Fig. 4 Heatmap illustrating the expression profile of differentially expressed lncRNAs (presented as Z-score values) in porcine granulosa cells treated with TCDD for 3, 12 or 24 h. Red blocks represent up- and green blocks represent down-regulated lncRNAs. The color scale of the heatmap shows the expression level where the brightest green stands for − 2.0 Z-score and the brightest red stands for + 2.0 Z-score. Z-score was calculated using logFPKM and the following equation: [x -mean(x)] / sd(x)where ‘sd’ means standard deviation. C 3–24: untreated porcine granulosa cells (control) cultured for3, 12or24 h.TCDD 3–24: porcine granulosa cells treated with TCDD (100 nmol/L) for 3, 12 or 24 h. *transcript expression was not detected in control (untreated) cells; ○ transcripts with expression level altered in two or three examined incubation times Among the identified events, 28 and 2 were defined as dif- porcine granulosa cells exposed to TCDD. The effects of ferentially expressed exons, whereas 5 and 2 were described TCDD on the porcine granulosa cell transcriptome have as splice junctions after 3 h and 12 h of TCDD treatment, been previously examined, however the studies have been respectively (see Additional file 7). The JunctionSet profile focused solely on the expression of genes (RNA-Seq [23]; plot for the exemplary lncRNA XLOC (XLOC_020835) is real-time PCR [29]). In comparison to humans and presented in Fig. 5. Two lncRNAs - TCONS_00040606 and rodents, information concerning the identification and TCONS_00040607 – maybeformedonthebasisof this expression of lncRNAs in pigs is limited [24, 50]. To the selected XLOC. Moreover, the expression of lncRNA best of our knowledge, lncRNAs have not yet been identi- TCONS_00040607 correlated negatively with the expres- fied in porcine granulosa cells and TCDD effects on the sion of AhR. lncRNA expression profile have also not been examined. A total of 1,666 lncRNAs and 42,913 mRNAs were Validation of RNA-Seq data by real-time PCR identified in porcine granulosa cells in the present study. To validate the RNA-Seq results, three up-regulated The majority of the identified lncRNAs were classified lncRNAs, i.e., TCONS_00016901, TCONS_00031035 and as intergenic lncRNAs. The identified lncRNAs were TCONS_00034713 were chosen for real-time PCR. The found to be shorter in comparison to those of mRNAs. expression of the selected lncRNAs entirely confirmed the In addition, the exon number of lncRNAs was lower and results obtained by RNA-Seq (Fig. 6). the exon length was longer than those of mRNAs. The features of the identified lncRNAs were similar to those Discussion found in other studies [24, 43, 51, 52]. Recent advances in sequencing technologies revealed The analysis of the expression profile of lncRNAs in that the genome is widely transcribed, yielding a large TCDD-treated porcine granulosa cells revealed the presence number of ncRNAs. In the current study, we investi- of 22 DELs. Similarly to previously reported TCDD-induced gated the genome-wide lncRNA expression profile of changes in granulosal gene expression [23], the most of Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 9 of 13 Fig. 5 Presentation of differential expression of exons (E) and the occurrence of splice junction (J) in the exemplary selected lncRNA XLOC (XLOC_020835) identified in porcine granulosa cells treated with TCDD. The upper panel shows the expression level estimates for the mean normalized read counts for each exon or splice junction of XLOC_020835 identified in TCDD-treated (blue) and control (red) cells. The lower panel shows the exonic regions (boxes, labelled E001-E003) and known splice junctions (solid line, labelled J004). Statistically significant differences (P-adjusted < 0.05) of usage exons and junctions are drawn in violet. C3: untreated porcine granulosa cells (control) cultured for 3 h. TCDD3: porcine granulosa cells treated with TCDD (100 nmol/L) for 3 h DELs (15 out of 22) were detected after 3 h of TCDD results of negatively co-expressed genes suggested a treatment. Since TCDD was shown to induce oxidative possible involvement of DELs in a variety of biological stress in the ovary [53, 54], the altered expression of processes such as “intracellular signal transduction”, lncRNAs after 3 h of TCDD treatment may be associ- “response to stimulus” and “negative regulation of bio- ated with the regulation of the early cellular stress re- logical process”. Interestingly, AhR was found among sponse to TCDD. these negatively trans-regulated genes. Although the Potential functions of DELs identified in porcine gran- downregulation of AhR is a typical response to TCDD, ulosa cells were indirectly predicted via searching their recent reports concerning ncRNAs may provide more target cis- and trans-regulated protein-coding genes. details on TCDD mechanism of action. Li et al. [55] Screening of genes which were located within 10 kb demonstrated that miR-203 (microRNA) suppressed the from DELs enabled us to select 10 cis-target genes, but expression of AhR in TCDD-treated human lung and none of them was enriched into GO classification. To hepatic cells. In the current study, we identified identify the trans-target genes presumptively regulated TCONS_00040607 i.e., lncRNA, the expression of which by lncRNAs, we correlated the expression levels of DELs correlated negatively with that of AhR after 3 h of TCDD and protein-coding genes. A total of 916 negative and treatment. It is possible that TCONS_00040607 affects the 1,143 positive correlations were detected. The GO expression of AhR and is involved in its negative Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 10 of 13 Fig. 6 Real-time PCR validation of three selected differentially expressed lncRNAs (TCONS_00016901, TCONS_00031035, TCONS_00034713) which were identified in TCDD-treated porcine granulosa cells (treated vs. untreated cells) by RNA-Seq. Data are expressed as mean ± SEM (n =3–4). Statistical analysis was performed using Student’s t-test. Asterisks designate statistical differences (P < 0.05). AU: arbitrary units; C: control; TCDD: 2,3,7,8-tetrachlorodibenzo-p-dioxin regulation during the cellular response to TCDD. These playing an important role in the metabolism of xenobi- results suggest that the regulation of AhR by lncRNAs otics [56], correlated positively with the expression of two may constitute part of cellular defense mechanism against DELs (TCONS_00034713 and TCONS_00031305). Due dioxins. This hypothesis, however, needs additional experi- to the fact that TCDD treatment usually induces CYP1A1 mental verification. expression, this enzyme is considered to be a molecular Some of the positively co-expressed genes were enriched marker of TCDD action [57, 58]. In our previous study, in GO terms associated with “cellular response to chem- the expression of CYP1A1 increased significantly in a ical stimulus” and “regulation of signal transduction”.We time-dependent manner after 3, 12 and 24 h of porcine found that the expression of CYP1A1, coding a protein granulosa cell incubation with TCDD [23]. We also Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 11 of 13 demonstrated that the TCDD binding to the porcine events of differential usage of exons and splice junctions in CYP1A1 active site resulted in a rapid closure of the en- lncRNAs loci were identified after 3 h and 12 h of cell zyme substrate channels. This phenomenon may partially incubation with TCDD. The majority of the identified explain TCDD’s high resistance to biodegradation [59]. If events were defined as differentially expressed exons the TCDD binding causes a continuous CYP1A1 blockage, and only a few events were described as splice junc- the cellular response to TCDD may induce an extended tions. It is of interest that two alternatively spliced synthesis of CYP1A1. The fact that the expression of two forms for XLOC_020835 were detected after TCDD DELs: TCONS_00034713 and TCONS_00031305 corre- treatment. In the current study, the expression of one lated positively with the expression of CYP1A1 indicates of these forms i.e., TCONS_00040607, was found to their supportive role in the cellular reaction to TCDD. negatively correlate with the expression of AhR.These Five DELs (TCONS_00005658; TCONS_00016901; facts implicate that TCONS_00040607 is involved in TCONS_00048979; TCONS_00060223; TCONS_0006 the regulation of the TCDD-affected AhR expression in 4401) were found to be expressed only in TCDD-treated porcine granulosa cells. cells. The expression of two DELs (TCONS_00048979 and TCONS_00060223) correlated negatively with the Conclusions expression of the same three genes: ERN1, LIMK2 and In the current study, we identified and characterized EPHB4. ERN1 is associated with endoplasmic reticulum lncRNAs of porcine granulosa cells. We also examined stress and ER protein folding [60, 61]. EPHB4, in turn, is the effects of TCDD on the lncRNA expression profile. linked with the proliferation of ovarian carcinoma cells The co-expression analysis revealed that the identified [62] and LIMK2 with actin microfilament disruption in lncRNAs may influence the expression of numerous porcine oocytes [63]. The obtained data imply that genes, including those involved in: 1) cellular response TCONS_00048979 and TCONS_00060223 are mediators to TCDD (e.g., AhR), 2) dioxin metabolism (e.g., of TCDD action in porcine granulosa cells. CYP1A1, MafF), 3) endoplasmic reticulum stress (e.g., In addition, we identified three DELs located in the ERN1) as well as 4) adhesion and proliferation of cells antisense strand of protein coding genes. These anti- (TGFβ1, NHS, LIMK2). Our results suggest that sense strands were found to positively or negatively lncRNAs constitute a part of the regulatory apparatus of regulate the expression of their sense counterparts and, TCDD action in porcine granulosa cells. This study may therefore, they attracted a lot of attention [64, 65]. In the provide the foundation for future research focused on current study, TCONS_00038918, TCONS_00030731 molecular effects exerted by TCDD in ovarian cells. and TCONS_00064964 were found to be located in the respective antisense strands of MafF, TGFβI and NHS. Additional files MafF belongs to the small MAF family of basic-leucin zipper transcription factors, which affect genes encoding Additional file 1: Summary of the mapping of the RNA sequences to the reference genome. (XLSX 13 kb) proteins responsible for xenobiotic metabolism and Additional file 2: Cis-target genes predicted to be regulated by antioxidation [66–68]. TGFβ1 encodes an extracellular differentially expressed lncRNAs in porcine granulosa cells exposed to matrix (ECM) protein reported to interact with various TCDD for 3, 12 and 24 h. (XLSX 50 kb) matrix molecules (collagens, fibronectin and laminin), Additional file 3: Negatively co-expressed trans-target genes predicted contributing to cell adhesion, proliferation and migration to be regulated by differentially expressed lncRNAs identified in porcine granulosa cells exposed to TCDD for 3, 12 and 24 h. (XLSX 78 kb) [69]. ECM proteins were found to be affected by TCDD Additional file 4: Positively co-expressed trans-target genes predicted in marmosets [70]. NHS products, in turn, were described to be regulated by differentially expressed lncRNAs identified in porcine to maintain cell morphology by remodeling the actin cyto- granulosa cells exposed to TCDD for 3, 12 and 24 h. (XLSX 126 kb) skeleton [71]. The in silico data concerning the possible Additional file 5: Functional enrichment analysis of the negatively relationships between lncRNAs and protein-coding genes co-expresed trans-target genes predicted to be regulated by differentially expressed lncRNAs identified in porcine granulosa cells exposed to TCDD in porcine granulosa cells treated with TCDD are sup- for 3, 12 and 24 h. (XLSX 54 kb) ported by the results of our recent study [72]. In this Additional file 6: Functional enrichment analysis of the positively study, the abundance of heat shock proteins as well as co-expresed trans-target genes predicted to be regulated by differentially cytoskeleton and ECM proteins were significantly affected expressed lncRNAs identified in porcine granulosa cells exposed to TCDD for 3, 12 and 24 h. (XLSX 36 kb) by TCDD in porcine granulosa cells. Additional file 7: Differential expression of exons and the occurence of Similarly to protein-coding genes, lncRNAs also undergo splice junction in XLOCs of the lncRNAs identified in porcine granulosa alternative splicing, resulting in the formation of numerous cells exposed to TCDD for 3 and 12 h. (XLSX 16 kb) lncRNA isoforms and, in consequence, extending their regulatory capabilities [73, 74]. It was demonstrated that Abbreviations xenobiotics may affect alternative splicing, modifying the AhR: Aryl hydrocarbon receptor; antisense lncRNAs: LncRNAs originate from process in a specific manner [75]. In the current study, the the opposite strand of DNA; ARNT: AhR nuclear translocator; bFGF: Fibroblast Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 12 of 13 growth factor-basic human; bHLH-PAS: The basic helix-loop-helix/PER-ARNT- 4. Liu D, Mewalal R, Hu R, Tuskan GA, Yang X. New technologies SIM; CNCI: Coding-Non-Coding Index; CPC: Coding Potential Calculator; accelerate the exploration of non-coding RNAs in horticultural plants. CYP1A1: Cytochrome P450 1A1; DELs: Differentially expressed lncRNAs; Hortic Res. 2017;4:17031. DMEM: Dulbecco’s modified Eagle’s medium; E : Estradiol; EPHB4: Ephrin 5. Prabhakar B, Zhong XB, Rasmussen TP. Exploiting long noncoding RNAs as type-B receptor 4; ERN1: Endoplasmic reticulum to nucleus signaling 1; pharmacological targets to modulate epigenetic diseases. Yale J Biol Med. FBS: Fetal bovine serum; FEELnc: FlExible Extraction of LncRNAs; GO: Gene 2017;90:73–86. Ontology database; LIMK 2: LIM domain kinase 2; lincRNAs: LncRNAs 6. Chen X, Yan CC, Zhang X, You ZH. Long non-coding RNAs and complex originate from intergenic regions of the genome; linRNAs: LncRNAs originate diseases: from experimental results to computational models. Brief from introns of protein-coding genes; lncRNAs: Long non-coding RNAs; Bioinform. 2016;18:558–76. MafF: MAF bZIP transcription factor F; ncRNAs: Non-coding RNAs; 7. Perry RB, Ulitsky I. The functions of long noncoding RNAs in development NEAA: MEM non-essential amino acid solution; NHS: NHS actin remodeling and stem cells. Development. 2016;143:3882–94. regulator; nt: Nucleotides; P : Progesterone; PAH: Polycyclic aromatic 4 8. Szcześniak MW, Makałowska I. lncRNA-RNA interactions across the human hydrocarbon; RIN: RNA integrity number; RNA-Seq: RNA sequencing; Transcriptome. PLoS One. 2016;11(3):e0150353. TCDD: 2,3,7,8-tetrachlorodibenzo-p-dioxin; TGFβI: Transforming growth factor 9. Guo Q, Cheng Y, Liang T, He Y, Ren C, Sun L, et al. Comprehensive beta-induced analysis of lncRNA-mRNA co-expression patterns identifies immune- associated lncRNA biomarkers in ovarian cancer malignant progression. Sci Rep. 2015;5:17683. Acknowledgements 10. Basham KJ, Leonard CJ, Kieffer C, Shelton DN, McDowell ME, Bhonde VR, We would like to thank Piotr Ciereszko for English language review. et al. Dioxin exposure blocks lactation through a direct effect on mammary epithelial cells mediated by the aryl hydrocarbon receptor repressor. Toxicol Funding Sci. 2015;14:36–45. This study was supported by National Science Centre (2012/05/B/NZ9/03333) 11. ATSDR. Draft Update Toxicological Profile for Chlorinated Dibenzo-p-dioxins. and The Ministry of Science and Higher Education in Poland (UWM No. Prepared by Research Triangle Institute for U.S. Department of Health and 528.0206.0806). HumanServices, Agency for Toxic Substances Disease Registry (ATSDR), Atlanta, GA. 677 pp +appendices (1998). Availability of data and materials 12. Nicolopoulou-Stamati P, Pitsos MA. The impact of endocrine disruptors on The datasets analyzed during the current study are available in the NCBI the female reproductive system. Hum Reprod. 2001;7:323–30. BioProject database under accession number: PRJNA429720 (https://www.n 13. Pocar P, Fischer B, Klonisch T, Hombach-Klonisch S. Molecular interactions of cbi.nlm.nih.gov/bioproject/?term=PRJNA429720). the aryl hydrocarbon receptor and its biological and toxicological relevance for reproduction. Reproduction. 2005;129:379–89. Authors’ contributions 14. Piekło R, Grochowalski A, Gregoraszczuk EL. 2,3,7,8-tetrachlorodibenzo- REC, AN, SS & LP contributed to the experimental concept; REC, MR, AN, SS, pdioxin alters follicular steroidogenesis in time- and cell-specific manner. AS & KO participated in creating the experimental design; MR, AS, AN & KO Exp Clin Endocr Diab. 2000;108:299–304. performed cell culture experiments and molecular experiments; LP & JPJ 15. Grochowalski A, Chrzaszcz R, Piekło R, Gregoraszczuk EL. Estrogenic and performed bioinformatics analysis of the data and participated in discussion antiestrogenic effect of in vitro treatment of follicular cells with 2,3,7,8- and interpretation of bioinformatics data; MR, AN & REC analyzed the data, tetrachlorodibenzo-p-dioxin. Chemosphere. 2001;43:823–7. discussed the obtained results and wrote manuscript; AS, KO, SS & TM 16. Petroff BK, Roby KF, Gao X, Son D, Williams S, Johnson D, et al. A participated in the interpretation of the data; MR, REC, AN & JPJ worked out review of mechanisms controlling ovulation with implications for the the figures and tables. All authors read and approved the manuscript. anovulatory effects of polychlorinated dibenzo-p-dioxin in rodents. Toxicology. 2001;158:91–107. 17. Moran FM, Vandevoort CA, Overstreet JW, Lasley BL, Conley AJ. Molecular Ethics approval and consent to participate target of endocrine disruption in human luteinizing granulosa cells by Not applicable. 2,3,7,8-tetrachlorodibenzo-p-dioxin: inhibition of estradiol secretion due to decrease 17a-hydroxylase/17,20-lyase cytochrome P450 expression. Consent for publication Endocrinology. 2003;144:467–73. Not applicable. 18. Mandal PK. Dioxin: a review of its environmental effects and its aryl hydrocarbon receptor biology. J Comp Physiol B. 2005;175:221–30. Competing interests 19. Jablonska O,Piasecka J,Petroff BK,Nynca A,Siawrys G, Wąsowska B, The authors declare that they have no competing interests. et al. In vitro effects of 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) on ovarian, pituitary, and pineal function in pigs. Theriogenology. 2011;76: Author details 921–32. Department of Animal Anatomy and Physiology, Faculty of Biology and 20. Albertini DF, Combelles CM, Benecchi E, Carabatsos MJ. Cellular basis for paracrine Biotechnology, University of Warmia and Mazury in Olsztyn, Oczapowskiego regulation of ovarian follicle development. Reproduction. 2001;121:647–53. 1A, 10-719 Olsztyn, Poland. Laboratory of Molecular Diagnostics, Faculty of 21. Jablonska O, Ciereszko RE. The expression of aryl hydrocarbon receptor in Biology and Biotechnology, University of Warmia and Mazury in Olsztyn, porcine ovarian cells. Reprod Domest Anim. 2013;48:710–6. Prawochenskiego 5, 10-720 Olsztyn, Poland. Department of Plant 22. Sadowska A, Nynca A, Korzeniewska M, Piasecka-Srader J, Jablonska M, Physiology, Genetics and Biotechnology, Faculty of Biology and Orlowska K, et al. Characterization of porcine Granulosa cell line AVG-16. Biotechnology, University of Warmia and Mazury in Olsztyn, Oczapowskiego Folia Biol (Praha). 2015;61:184–94. 1A, 10-719 Olsztyn, Poland. 23. Sadowska A, Nynca A, Ruszkowska M, Paukszto L, Myszczynski K, Orlowska K, et al. Transcriptional profiling of porcine granulosa cells exposed to 2,3,7,8- Received: 7 May 2018 Accepted: 23 August 2018 tetrachlorodibenzo-p-dioxin. Chemosphere. 2017;178:368–77. 24. Tang Z, Wu Y, Yang Y, Yang YT, Wang Z, Yuan J, et al. Comprehensive analysis of long non-coding RNAs highlights their spatio-temporal expression patterns and evolutional conservation in Sus scrofa. Sci Rep. References 2017;7:43166. 1. Dempsey JL, Cui JY. Long non-coding RNAs: a novel paradigm for toxicology. Toxicol Sci. 2017;155:3–21. 25. Horisberger MA. A method for prolonged survival of primary cell lines. In 2. Yu L, Tai L, Zhang L, Chu Y, Li Y, Zhou L. Comparative analyses of long non- Vitro Cell Dev Biol Anim. 2006;42:143–8. coding RNA in lean and obese pigs. Oncotarget. 2017;8:41440–50. 26. Gregoraszczuk EL, Wojtowicz AK, Zabielny E, Grochowalski A. Dose-and time 3. Dykes IM, Emanueli C. Transcriptional and post-transcriptional gene dependent effect of 2,3,7,8-tetrachlorodibenzo-p-dioxin (TCDD) on regulation by long non-coding RNA. Genomics Proteomics Bioinformatics. progesterone secretion by porcine luteal cells cultured in vitro. J Physiol 2017;5:177–86. Pharmacol. 2000;51:127–35. Ruszkowska et al. Journal of Animal Science and Biotechnology (2018) 9:72 Page 13 of 13 27. Gregoraszczuk EL. Dioxin exposure and porcine reproductive hormonal 52. Xia J, Xin L, Zhu W, Li L, Li C, Wang Y, et al. Characterization of long activity. Cad Saude Publ. 2002;18:453–62. noncoding RNA transcriptome in highenergy diet induced non-alcoholic 28. Jablonska O, Piasecka-Srader J, Nynca A, Kołomycka A, Robak A, Wąsowska steatohepatitis minipigs. Sci Rep. 2016;6:30709. B, et al. 2,3,7,8-Tetrachlorodibenzo-p-dioxin alters steroid secretion but does 53. Melekoglu R, Ciftci O, Cetin A, Basak N, Celik E. The beneficial effects of not affect cell viability and the incidence of apoptosis in porcine luteinised Montelukast against 2,3,7,8-tetrachlorodibenzo-p-dioxin toxicity in female granulosa cells. Acta Vet Hung. 2014;62:408–21. reproductive system in rats. Acta Cir Bras. 2016;31:557–63. 29. Piasecka-Srader J, Sadowska A, Nynca A, Orlowska K, Jablonska M, Jablonska 54. Bhattacharya P, Keating AF. Impact of environmental exposures on ovarian O, et al. The combined effects of 2,3,7,8-tetrachlorodibenzo-p-dioxin and function and role of xenobiotic metabolism during ovotoxicity. Toxicol Appl the phytoestrogen genistein on steroid hormone secretion, AhR and ERβ Pharmacol. 2017;261:227–35. expression and the incidence of apoptosis in granulosa cells of medium 55. Li D, Liu C, Yu H, Zeng X, Xing X, Chenet L, et al. AhR is negatively porcine follicles. J Reprod Dev. 2016;62:103–13. regulated by miR-203 in response to TCDD or BaP treatment. Toxicol Res. 2014;3:142–51. 30. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina 56. Androutsopoulos VP, Tsatsakis AM, Spandidos DA. Cytochrome P450 CYP1A1: sequence data. Bioinformatics. 2014;30:2114–20. wider roles in cancer progression and prevention. BMC Cancer. 2009;9:187. 31. Dobin A, Gingeras TR. Mapping RNA-seq reads with STAR. Curr Protoc 57. Whitlock JP Jr. Induction of cytochrome P4501A1. Annu Rev Pharmacol Bioinformatics. 2015;51:11.14.1–19. Toxicol. 1999;39:103–25. 32. Trapnell C, Roberts A, Goff L, Pertea G, Kim D, Kelley DR, et al. Differential 58. Kim S, Dere E, Burgoon LD, Chang CC, Zacharewski TR. Comparative analysis gene and transcript expression analysis of RNA-seq experiments with of AhR-mediated TCDD-elicited gene expression in human liver adult stem TopHat and cufflinks. Nat Protoc. 2012;7:562–78. cells. Toxicol Sci. 2009;112:229–44. 33. Pertea M, Pertea GM, Antonescu CM, Chang TC, Mendell JT, Salzberg SL. 59. Molcan T, Swigonska S, Orlowska K, Myszczynski K, Nynca A, Sadowska A, StringTie enables improved reconstruction of a transcriptome from RNA-seq et al. Structural-functional adaptations of porcine CYP1A1 to metabolize reads. Nat Biotechnol. 2015;33:290–5. polychlorinated dibenzo-p-dioxins. Chemosphere. 2017;168:205–16. 34. Pertea M, Kim D, Pertea GM, Leek JT, Salzberg SL. Transcript-level expression 60. Minchenko OH. Inhibition of ERN1 signalling enzyme affects hypoxic analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat regulation of the expression of E2F8, EPAS1, HOXC6, ATF3, TBX3 AND Protoc. 2016;11:1650–67. FOXF1 genes in U87 glioma cells. Ukr Biochem J. 2015;87:76–87. 35. Zhan S, Dong Y, Zhao W, Guo J, Zhong T, Wang L, et al. Genome-wide 61. Rashid HO, Yadav RK, Kim HR, Chae HJ. ER stress: autophagy induction, identification and characterization of long non-coding RNAs in inhibition and selection. Autophagy. 2015;11:1956–77. developmental skeletal muscle of fetal goat. BMC Genomics. 2016;7:666. 62. Kumar SR, Masood R, Spannuth WA, Singh J, Scehnet J, Kleiber G, et al. The 36. Han DX, Sun XL, Fu Y, Wang CJ, Liu JB, Jiang H, et al. Identification of long receptor tyrosine kinase EphB4 is overexpressed in ovarian cancer, provides non-coding RNAs in the immature and mature rat anterior pituitary. Sci Rep. survival signals and predicts poor outcome. Br J Cancer. 2007;96:1083–91. 2017;7:17780. 63. Jia RX, Duan X, Song SI, Su SC. LIMK1/2 inhibitor LIMKi 3 suppresses porcine 37. Kong L, Zhang Y, Ye ZQ, Liu XQ, Zhao SQ, Wei L, et al. CPC: assess the oocyte maturation. PeerJ. 2016;4:e2553. protein-coding potential of transcripts using sequence features and support 64. Khorkova O, Myers AJ, Hsiao J, Wahlestedt C. Natural antisense transcripts. vector machine. Nucleic Acids Res. 2007;35:W345–9. Hum Mol Genet. 2014;23(R1):R54–63. 38. Sun L, Luo H, Bu D, Zhao G, Yu K, Zhang C, et al. Utilizing sequence intrinsic 65. Villegas VE, Zaphiropoulos PG. Neighboring gene regulation by antisense composition to classify protein-coding and long non-coding transcripts. long non-coding RNAs. Int J Mol Sci. 2015;16:3251–66. Nucleic Acids Res. 2013;41(17):e166. 66. Dhakshinamoorthy S, Jaiswal AK. Small maf (MafG and MafK) proteins 39. Wucher V, Legeai F, Hédan B, Rizk G, Lagoutte L, Leeb T, et al. FEELnc: a tool negatively regulate antioxidant response element-mediated expression and for long non-coding RNA annotation and its application to the dog antioxidant induction of the NAD(P)H:Quinone oxidoreductase1 gene. J Biol transcriptome. Nucleic Acids Res. 2017;45(8):e57. Chem. 2000;275:40134–41. 40. Mistry J, Bateman A, Finn RD. Predicting active site residue annotations in 67. Katsuoka F, Motohashi H, Ishii T, Aburatani H, Engel JD, Yamamoto M. Genetic the Pfam database. BMC Bioinformatics. 2007;8:298. evidence that small maf proteins are essential for the activation of antioxidant 41. Li A, Zhang J, Zhou Z. PLEK: a tool for predicting long non-coding RNAs response element-dependent genes. Mol Cell Biol. 2005;25:8044–51. and messenger RNAs based on an improved k-mer scheme. BMC 68. Massrieh W, Derjuga A, Doualla-Bell F, Ku CY, Sanborn BM, Blank V. Bioinformatics. 2014;15:311. Regulation of the MAFF transcription factor by proinflammatory cytokines in 42. Kung JT, Colognori D, Lee JT. Long noncoding RNAs: past, present, and myometrial cells. Biol Reprod. 2006;74:699–705. future. Genetics. 2013;193:651–69. 69. Kim JE. Molecular properties of wild-type and mutant betaIG-H3 proteins. 43. Huber W, Carey VJ, Gentleman R, Anders S, Carlson M, Carvalho BS, et al. Invest Ophthalmol Vis Sci. 2002;43:656–61. Orchestrating high-throughput genomic analysis with bioconductor. Nat 70. Riecke K, Grimm D, Shakibaei M, Kossmehl P, Schulze-Tanzil G, Paul M, et al. Methods. 2015;12:115–21. Low doses of 2,3,7,8-tetrachlorodibenzo-p-dioxin increase transforming 44. Ren H, Wang G, Chen L, Jiang J, Liu L, Li N, et al. Genome-wide analysis of growth factor beta and cause myocardial fibrosis in marmosets long non-coding RNAs at early stage of skin pigmentation in goats (Capra (Callithrixjacchus). Arch Toxicol. 2002;76:360–6. hircus). BMC Genomics. 2016;17:67. 71. Zhang DD, Du JZ, Topolewski J, Wang XM. Review recent progress in 45. Zhang T, Zhang X, Han K, Zhang G, Wang J, Xie K, et al. Genome-wide identification and characterization of loci associated with sex-linked analysis of lncRNA and mRNA expression during differentiation of congenital cataract. Genet Mol Res. 2016;29, 15(3) abdominal Preadipocytes in the chicken. G3 (Bethesda). 2017;7:953–66. 72. Orlowska K, Swigonska S, Sadowska A, Ruszkowska M, Nynca A, Molcan T, 46. Reimand J, Kull M, Peterson H, Hansen J, Vilo J. G:profiler—a web-based et al. The effects of 2,3,7,8-tetrachlorodibenzo-p-dioxin on the proteome of toolset for functional profiling of gene lists from large-scale experiments. porcine granulosa cells. Chemosphere. 2018;212:170–81. Nucleic Acids Res. 2007;35:W193–200. 73. Bozgeyik E, Igci YZ, Sami Jacksi MF, Arman K, Gurses SA, Bozgeyik I, et al. A 47. Hartley SW, Mullikin JC. QoRTs: a comprehensive toolset for quality novel variable exonic region and differential expression of LINC00663 non- control and data processing of RNA-Seq experiments. BMC Bioinformatics. coding RNA in various cancer cell lines and normal human tissue samples. 2015;16:224. Tumour Biol. 2016;37:8791–8. 48. Hartley SW, Mullikin JC. Detection and visualization of differential splicing in 74. Sen R, Doose G, Stadler PT. Rare splice variants in long non-coding RNAs. RNA-Seq data with JunctionSeq. Nucleic Acids Res. 2016;44:e127. Non-coding RNA. 2017;3(3):23. 49. Hellemans J,Mortier G, De PaepeA,Speleman F,VandesompeleJ. qBase relative 75. Zaharieva E, Chipman JK, Soller M. Alternative splicing interference by quantification framework and software for management and automated analysis xenobiotics. Toxicology. 2012;296:1–12. of real-time quantitative PCR data. Genome Biol 2007;8(2):R19. 50. Zhao W, Mu Y, Ma L, Wang C, Tang Z, Yang S, et al. Systematic identification and characterization of long intergenic non-coding RNAs in fetal porcine skeletal muscle development. Sci Rep. 2015;5:8957. 51. Ran M, Chen B, Li Z, Wu M, Liu X, He C, et al. Systematic identification of long noncoding RNAs in immature and mature porcine testes. Biol Reprod. 2016;94:77.

Journal

Journal of Animal Science and Biotechnology – Springer Journals

Published: Oct 11, 2018

There are no references for this article.