Transcriptomic Changes in Response to Form of Selenium on the Interferon-Tau Signaling Mechanism in the Caruncular Tissue of Beef Heifers at Maternal Recognition of Pregnancy

We have reported that selenium (Se) provided to grazing beef cattle in an inorganic (ISe) form versus a 1:1 mixture (MIX) of inorganic and organic (OSe) forms affects cholesterol biosynthesis in the corpus luteum (CL), the abundance of interferon tau (IFNτ) and progesterone (P4)-induced mRNAs in the caruncular (CAR) tissue of the endometrium, and conceptus length at maternal recognition of pregnancy (MRP). In this study, beef heifers were supplemented with a vitamin–mineral mix containing 35 ppm Se as ISe or MIX to achieve a Se-adequate status. Inseminated heifers were killed at MRP (d 17, n = 6 per treatment) for tissue collection. In CAR samples from MIX versus ISe heifers, qPCR revealed that mRNA encoding the thyroid regulating DIO2 and DIO3 was decreased (p < 0.05) and a complete transcriptomic analysis revealed effects on the interferon JAK-STAT1/2 pathway, including decreased expression of mRNAs encoding the classical interferon stimulated genes IFIT1, IFIT2, IFIT3, IRF1, IRF9, ISG15, OAS2, and RSAD2 (p < 0.05). Treatment also affected the abundance of mRNAs contributing to the immunotolerant environment (p < 0.05). In combination, these findings suggest more advanced preparation of the CAR and developing conceptus for implantation and to evade immune rejection by the maternal system in MIX- vs. ISe-treated heifers.

Across the United States, the concentration of Se in the soils and, subsequently, in the forages available to grazing livestock varies by geographic region [19].The concentration of Se in the forages and grains in the southeastern United States, including Kentucky, are low, at less than 0.05 ppm, to variable.Only about 50% of the feedstuffs produced contain greater than 0.1 ppm [19], the minimal recommended dietary intake of this trace mineral [20].This provides a challenge for cattle producers in these Se-deficient regions; hence, it is necessary to supplement Se in the diet of cattle grazing these forages [21][22][23][24][25]. Providing this auxiliary Se to overcome a deficiency in this trace mineral is essential for various physiological processes, including immune function [24,26], growth [22], antioxidant status [27], and fertility [23].Commonly, supplemental Se is provided in a vitamin-mineral mix as an inorganic form (ISe), sodium selenite, or sodium selenate.However, organic forms of selenium (OSe), and particularly selenomethionine, are available when cattle graze forage [19].
We consistently studied the effects of form of supplemental Se as either ISe or a 50%:50% mixture of ISe to OSe (MIX) on gonadal function in cattle and reported a MIXinduced increase in the systemic concentration of progesterone (P4) in the early luteal phase of the estrous cycle [28,29], with this increase occurring at a time that can significantly improve embryonic development by altering endometrial function that supports conceptus elongation [30][31][32] and the production of interferon tau (IFNτ), a protein from the developing conceptus that signals maternal recognition of pregnancy (MRP [33]).Indeed, we observed a significant increase in the length of the conceptus in MIX versus ISe-supplemented heifers, which has clear implications for fertility outcomes as it may be better prepared for the impending processes of attachment and implantation [34].Receptivity of the endometrium to allow implantation of the conceptus is dependent on P4 [35], uterine restructuring to promote placentation [36], and the ability of the fetal allograft to evade the maternal immune response [37].Considerably, we performed a targeted qPCR analysis, determining changes in the relative abundance of several P4 and IFNτregulated mRNA transcripts in the caruncular (CAR) tissues of the endometrium from ISeversus MIX-supplemented heifers at MRP (d 17 of gestation) and observed a significant decrease in the expression of mRNAs encoding diacylglycerol o-acyltransferase 2 (DGAT2), fibroblast growth factor 2 (FGF2), interferon induced protein with tetratricopeptide repeats 3 (IFIT3), ISG15 ubiquitin like modifier (ISG15), MX dynamin like GTPase 1 (MX1), 2 -5oligoadenylate synthetase 2 (OAS2), and radical S-adenosyl methionine domain containing 2 (RSAD2 [34]), with this study expanding on those results.
The necessary immunotolerant environment is the result of a careful balance between compositional changes in the innate immune system, which appears to remain active to protect against pathogens [38], and the adaptive immune system, which must be repressed [37,38].During this time period, ligands from the conceptus signal the endometrium to effect necessary changes, and in the presence of P4, there is an increase in macrophages, dendritic cells, and NK cells that migrate into uterine epithelial, endothelial, and stomal cells [37], as well as changes in the abundance of molecules that promote immune tolerance such as indoleamine 2,3-dioxygenase 1 (IDO1), transforming growth factor beta (TGFB), and interferon stimulated genes (ISG) [37].
In our present experiment, commercial Angus-cross heifers (N = 20) received 45-day ad libitum access to a basal vitamin-mineral mix with no exogenous source of Se (depletion phase), followed by 45-day supplementation with a vitamin-mineral mix containing 35 ppm Se as ISe (repletion phase), and then were randomly assigned to received ad libitum access to supplemental treatment as either inorganic (ISe; sodium selenite, n = 10) or a 1:1 combination (MIX, n = 10) of ISe and OSe (1:1 sodium selenite and SEL-PLEX) for at least 90 days prior to estrous synchronization and artificial insemination.Heifers were then killed at d 17 of presumed pregnancy.An intact conceptus and other experimental samples were collected from both ISe-and MIX-treated heifers (n = 6 per treatment) immediately after death.Caruncular (CAR) tissue was used herein; however, the embryo, CAR and intercaruncular (ICAR) tissues, whole blood, and serum were utilized in Crites et al. [34].
This experiment had two primary objectives.Firstly, we aimed to quantify the changes in mRNA transcripts of selenoproteins with known or proposed antioxidant capabilities and select selenoprotein receptors in response to the form of Se at MRP.Our hypothesis was that a cohort of these enzymes, and particularly the cytosolic GPX1, will be of greater abundance in MIX-vs.ISe-treated heifers, which would correspond with the whole animal's ability to mitigate the increased ROS that is associated with an increased metabolic demand at MRP. Subsequently, with the absence of a global perspective on the effect of form of Se on mRNA transcripts in the caruncular regions of the endometrium, the second objective of this study was to examine transcriptomic changes in the CAR tissue at MRP, postulating that we would observe MIX-induced changes in transcript expression indicating advancement of the CAR transcriptome in preparation for the impending processes of attachment and implantation.Defining the mRNA profile of selenoproteins and the transcriptome of CAR tissue as a whole will provide insight into how the form of Se is related to improved fertility and longer conceptuses at MRP.Ultimately, this provides a producer-friendly and sustainable method to increase whole farm productivity and profitability of grazing beef cattle.

Real-Time PCR Analysis of Selenoproteins and Selenoprotein P Receptors mRNA
A total of twenty-nine selenoprotein and selenoprotein P receptor transcripts were targeted via qPCR analysis.The MIX form of supplemental Se significantly (p < 0.05) decreased the relative abundance of Dio2 and Dio3 mRNAs (Table 1).The relative abundance of transcripts encoding the other targeted selenoproteins and the selenoprotein P receptors was not affected by form of Se.RNAseq also indicated that the MIX form of Se significantly (p < 0.05) decreased the relative abundance of Txnrd2 and Selenon mRNAs, and increased the abundance of Selnok mRNA; however, these findings were not confirmed by qPCR.

Cluster Analyses
Principal component analysis (PCA) of all RNA-Seq data was performed to evaluate the relative relationships and variation among the individual heifers.The score plot (Figure 1A) demonstrates that principal component 1 (PC#1, x-axis) explained 25% of variance among the samples and principal component 2 (PC#2, y-axis) explained 12% of variance.Results indicate that there is slight separation between treatment groups, with some similarities in the expression of DEGs.

Cluster Analyses
Principal component analysis (PCA) of all RNA-Seq data was performed to the relative relationships and variation among the individual heifers.The score p ure 1A) demonstrates that principal component 1 (PC#1, x-axis) explained 25% of among the samples and principal component 2 (PC#2, y-axis) explained 12% of v Results indicate that there is slight separation between treatment groups, with so ilarities in the expression of DEGs.
Hierarchical clustering analysis of the top 1000 DEGs (Figure 1B) indicates a erable separation between Se treatments; however, there may be some overlap scriptomic profiles of the CAR from ISe-and MIX-supplemented heifers.Hierarchical clustering analysis of the top 1000 DEGs (Figure 1B) indicates a considerable separation between Se treatments; however, there may be some overlap in transcriptomic profiles of the CAR from ISe-and MIX-supplemented heifers.

Differentially Expressed Genes
The Wald test of DESeq2 was used to identify changes in CAR between the ISe-and MIX-supplemented heifers.At p < 0.05, there were 2029 DEGs with 1038 transcripts upregulated, and 991 transcripts downregulated in MIX compared to ISe.The most differentiated genes that were increased and decreased according to fold change in MIX compared to ISe are provided in Table 2.

Pathway and Gene Network Analysis
To determine global changes in DEGs in response to form of Se, bioinformatic analysis was performed using QIAGEN's Ingenuity Pathway Analysis (IPA, QIAGEN, Redwood City, CA, USA).Canonical pathway analysis indicated that the top five canonical pathways (Table 3) based on p-values affected by the form of Se are ATM signaling (p < 0.0001), spliceosomal cycle (p < 0.0001), role of protein kinase R (PKR) in interferon induction and antiviral response (p < 0.0001), coordinated lysosomal expression and regulation (CLEAR) signaling pathway (p < 0.0001), and autophagy (p < 0.0001).Further, Figure 2 shows the top canonical pathways in CAR ranked by z-score, which is used to predict activation state.Interestingly, the most significantly affected pathways that were positively affected are spliceosomal cycle (z-score = 2.673), cold shock domain containing E1 (CSDE1) signaling pathway (z-score = 2.138), cell cycle control of chromosomal replication (z-score = 1.508), cell cycle: G1/S checkpoint regulation (z-score = 1.414), and 3-phosphoinositide biosynthesis (z-score = 1.342).The most negatively affected pathways based on z-score are the role of PKR in interferon induction and antiviral response (z-score = −3.128),interferon signaling (z-score = −2.333),CD40 signaling (z-score = −2.333),death receptor signaling (z-score = −2.183),and TNF receptor superfamily member 1A (TNFR1) signaling (z-score = −2.121).
The top upstream regulators identified using IPA at p < 0.05 were hepatocyte nuclear factor 4 alpha (HNF4A), estrogen receptor 1 (ESR1), KRAS proto-oncogene, GTPase (KRAS), transforming growth factor beta 1 (TGFB1), and Pyridostatin.The top five molecular and cellular functions with specific actions identified by IPA are indicated in Table 4.Only functions with a z-score greater than or equal to the absolute value of 0.5 are reported. 1Selenium was supplemented at 35 ppm as either inorganic (ISe; sodium selenite) or a 1:1 combination (MIX) of ISe and OSe (1:1 sodium selenite and SEL-PLEX).Selenium was supplemented to treatment groups ad libitum. 2 Results from IPA were obtained based on log2 fold changes calculated using DESeq2. 3The ratio calculated as the number of differentially expressed genes (p < 0.05) in a given pathway divided by the total number of genes that make up that pathway. 4p-values were determined using IPA.
Formation of cellular protrusions −0.647 <0.0001 Replication of centriole 0.921 <0.0001 1 Selenium was supplemented at 35 ppm as either inorganic (ISe; sodium selenite) or a 1:1 combination (MIX) of ISe and OSe (1:1 sodium selenite and SEL-PLEX).Selenium was supplemented to treatment groups ad libitum. 2 Results from IPA were obtained based on log2 fold changes calculated using DESeq2. 3Z-score was calculated using IPA and was used to infer activation state.Only fractions resulting in a z-score with an absolute value of ≥0.5 are reported herein. 4p-values were determined using IPA.   1 Selenium was supplemented at 35 ppm as either inorganic (ISe; sodium selenite) or a 1:1 combination (MIX) of ISe and OSe (1:1 sodium selenite and SEL-PLEX).Selenium was supplemented to treatment groups ad libitum. 2 Results from IPA were obtained based on log2 fold changes calculated using DESeq2. 3Z-score was calculated using IPA and was used to infer activation state.Only fractions resulting in a z-score with an absolute value of ≥0.5 are reported herein. 4p-values were determined using IPA.

Real-Time PCR Analysis of Select mRNA Transcripts
To validate the outcomes of the RNA-sequencing, real-time PCR analysis was conducted on select transcripts that are responsive to interferons or are relative to uterine function at MRP (Table 5).The results of qPCR were consistent with RNA-sequencing analysis at either significance, declared at p ≤ 0.05, or a tendency at 0.05 < p ≤ 0.10, except for Irf9, Isg20, and Scarb1 mRNA.However, these were consistent in trend when compared to RNA-sequencing results, but were not significant for qPCR analyses.We did not corroborate all identified RNA-sequencing transcripts of interest; some are included below for purposes of clarity in the discussion.

Discussion
Soils and hence forages in geographic regions that are low in Se require producers to supplement this trace mineral in the diet of beef cattle to overcome Se-deficiencyassociated challenges with immune function [24,26], growth [22], antioxidant status [27], and fertility [23].Commonly, supplemental Se is provided in a vitamin-mineral mix as ISe; however, Ose forms are available when cattle graze forage [19].
The organic forms of Se are more bioavailable than the inorganic forms [39], and organic selenium can raise whole blood levels of Se more effectively compared to the inorganic sodium selenite or sodium selenate [40].Additionally, transcriptomic analysis in the liver of growing beef heifers supplemented with Se as either ISe, OSe, or a 50%:50% mixture of ISe or OSe (MIX) demonstrated that the form of Se creates specific phenotypes respective to form, with distinct physiological capabilities including changes in redox potential.Additionally, these transcriptomic changes were not indicative of a gradient effect respective to the different forms of Se utilized [41].
Using an extensively reported model [10,28,29,42,43], we studied the effects of supplemental form of Se as ISe (the industry standard) or a 1:1 mixture of ISe to OSe (MIX) on grazing beef cows and heifers and observed significant effects on the CL [10] and endometrium [34], which laid the premise for the present study.This experiment had two specific objectives.Firstly, we aimed to quantify the changes in mRNA transcripts of selenoproteins with known or proposed antioxidant capabilities and select selenoprotein P receptors in caruncular (CAR) tissue in response to form of Se at MRP.Secondly, this study sought to augment previous findings in our lab that identified key caruncular changes in the abundance of mRNA encoding several P4 and IFNτ-induced endometrial transcripts, which occurred synonymously with longer conceptuses in MIX-versus ISe-supplemented heifers at MRP [34].Therefore, we performed a transcriptomic analysis to elucidate global changes in the CAR transcriptome in response to form of Se.These results expand our understanding of mechanistic changes in the CAR that may be contributing to the previously observed significant differences in conceptus length at MRP and may be better preparing the CAR for attachment and implantation.

Selenoproteins
There are 25 selenoprotein genes identified in mammals [8,9], and we sought to analyze the Se-form effects on mRNA transcripts, encoding these, as well as select selenoprotein P receptors, in the caruncular tissue of heifers at MRP (d 17 of gestation).We hypothesized that a cohort of these enzymes, and particularly the cytosolic GPX1, will be of greater abundance in MIX-vs.ISe-treated heifers, which would indicate changes in the ability to mitigate potentially elevated reactive oxygen species associated with an increased metabolic demand at MRP, with the CAR at this time differentiating in early preparation for implantation and placentation.
Unexpectedly, we only observed changes in the abundance of mRNA transcripts encoding the intracellular iodothyronine deiodinases 2 and 3 (Dio2 and Dio3), and the relative abundance of mRNA encoding these two selenoproteins was significantly lower (p < 0.05) in MIX-supplemented heifers compared to those on ISe alone.Iodothyronine deiodinases (DIOs) are a family of selenoproteins that regulate the activation and deactivation of thyroid hormones, and, subsequently, growth, development, thermogenesis, cell differentiation and proliferation, and energy metabolism, including regulating metabolism of carbohydrates, proteins, and lipids [44,45].Regulating the circulating and intracellular concentrations of thyroid hormones is critical for placental development and function, as well as fetal growth [45].
DIO1 is bound in the plasma membrane and has been primarily localized to the liver, kidney, and thyroid [8].This selenoprotein converts T4 to T3 by deiodination of the outer ring, and it can convert T3 and T4 into the inactive T2 or rT3, respectively, via inner ring deiodination [46].Given the primary localization and function of DIO1, it is not unsurprising that we were unable to detect a measurable concentration of mRNA encoding this protein in CAR samples; however, it should be noted that previous studies have shown DIO1 to be abundant in the term placenta [47], which may be associated with circulating concentrations of T3 being highest at parturition [47,48].
Intracellular DIO2 is located in the endoplasmic reticulum and is primarily expressed in the nervous system, pituitary, thyroid, and brown adipose tissue [8].Low levels of serum T4 can result in an elevation in the activity of DIO2 [49].DIO3 is another membranebound selenoprotein that is primarily expressed in vascular tissue, skin, and the placenta, as well as being highly expressed in fetal and neonatal tissues [44].Its high level of expression in the placenta is to protect the growing fetus from the activity of maternal thyroid hormone as DIO3 primarily inactivates T3 [50].However, in ruminants, the placenta may be impermeable to transfer of thyroid hormones [45].It seems plausible that T3 could be effectively regulating the structural and metabolic changes occurring in the CAR in preparation for placental development in our experimental animals, which would be consistent with the previously observed longer embryos in MIX-versus ISe-supplemented heifers at MRP [34].

Caruncular Transcriptomics
With the absence of a global perspective of mRNA transcripts in the CAR in response to form of Se, we performed a transcriptomic analysis using next-generation sequencing to reveal the form of Se-affected pathways that may advance or impede the impending placental attachment.Our previous observation that MIX-supplemented heifers had significantly longer conceptuses on d 17 of gestation, coincident with a decreased level of expression of known P4-and INF-induced transcripts, including DGAT2, MX1, and OAS2 in the CAR [34], suggests that the process of MRP is more advanced in these heifers, as may be the development of the caruncular transcriptome in preparation for the initial process of adhesion and attachment, which is observed as early as days 16-18 of pregnancy in cattle [51,52].Results from the current transcriptomic analysis correspond with this, as evidenced by the significant changes in both interferon signaling, which is affecting the transcription of interferon responsive genes, and by the noted changes in transcripts regulating cellular organization, survival, and death.
Prominently, the top two canonical pathways downregulated in MIX-treated heifers were associated with endometrial response to type I interferon signaling ("Role of PKR in interferon induction and antiviral response" and "Interferon signaling").The trophoblastderived IFNτ is a type I interferon that signals by binding to the ubiquitously expressed type 1 interferon receptors (IFNAR1 and IFNAR2) on the endometrial luminal and superficial glandular epithelium [38,53].The primary function of this signal is to inhibit transcription of the estrogen (ER), effectively preventing estrogen-induced increases in the synthesis of ER, nuclear progesterone receptor (PGR), and oxytocin receptor (OTR), thus blocking production of luteolytic pulses of prostaglandin F 2α (PGF 2α ).This mechanism also appears to involve the actions of IRF2 [30].Presently, and reported in Crites et al. [34], the presence of a more developed conceptus in MIX-versus ISe-treated heifers demonstrated successful signaling at MRP.However, the absence of Se-form effects on the abundance of mRNA encoding IFNAR1, IFNAR2, IRF2, ER, or OXTR suggests that there are mechanisms still to be elucidated that may be of significance.Additionally, as reported in Crites et al. [34], qPCR analysis revealed a tendency for Pgr mRNA to be less abundant in MIX, which is consistent with the previously observed decrease in the expression of Pgr in response to continuous exposure of P4 in the endometrium [54].
In conjunction with being antiluteolytic, IFNτ exerts immense effects in the endometrium during MRP by promoting the transcription of interferon stimulated genes (ISGs, [55,56]).To do so, IFNτ binds IFNAR1 and IFNAR2 to activate the JAK/STAT signal transduction pathway.STAT1 can form a homodimer that translocates into the nucleus to stimulate the transcription of specific genes by binding to interferon-gamma-activated sequence (GAS).Alternatively, STAT1 and STAT2 form a heterodimer and interact with IRF9 to form the transcriptionally active gene complex factor 3. This complex migrates into the nucleus to bind interferon stimulated response elements (ISRE) that facilitate transcription of ISGs (Figure 3) [38,56,57].
conceptus occurred earlier in MIX supplemented heifers, as there was a stark decrease i the abundance of mRNA transcripts encoding IRF9.Notably, Irf9 mRNA is typically mor abundant between days 14-18 compared to days 25-40 of pregnancy in cattle, coinciden with decreases in mRNA encoding IRF3, MX1α, MX1β, and MX2 in CAR at these two tim points [59].To corroborate these findings, quantifying IFNτ in the uterine fluid would provide an indication of the concentration of this protein in the histotroph.Unfortunately in the present study, recovery of the uterine fluid at sample collection proved too incon sistent for subsequent quantification of the concentrations of IFNτ protein.For successful attachment and implantation to occur, there must also be restructurin of the endometrium and evasion of the maternal immune system by the developing con ceptus [36,37].Functional analysis of the significantly affected mRNA transcripts in th present study identified significant changes in the "spliceosomal cycle", "cell death and survival", "cellular assembly and organization", and "autophagy" pathways, suggestin Mechanistic illustration of downregulation of interferon-responsive genes and JAK/STAT1/2 signaling pathway.Red arrows demonstrate lower transcript abundance in MIX-Secompared to ISe-supplemented heifers.All transcripts included were significant via transcriptomic analysis, and an asterisk (*) indicates those transcripts that were corroborated via qPCR at p ≤ 0.05.Results of qPCR analysis for MX1 are reported in Crites et al. [34].
Results of the present transcriptomic analysis revealed significant MIX-form downregulation in the abundance of ISGs, particularly to those that are activated by ISRE.During RNA-sequencing analysis, we observed a significant MIX-induced decrease in Stat2 and Irf9 mRNA with no change in mRNA encoding STAT1.This suggests that downregulation of ISG is occurring via the interferon gene complex factor 3 protein complex that translocates into the nucleus to bind ISRE in the promoter regions of ISGs [38].Additionally, STAT2 and IRF9 self-regulate transcription of their own genes, thus leading to potentially further effects on the type 1 interferon JAK/STAT signaling transduction pathway [58].
Of the classical ISGs that are transcribed following activation of ISRE, mRNA transcripts for IFIT2, IFIT3, IRF9, ISG15, and STAT2 were significantly downregulated in MIX-supplemented heifers (Figure 3).We further detected significant decreases in mRNAs encoding the interferon-stimulated genes OAS2, RSAD2 and ISG20, and we previously reported a MIX-Se downregulation of Mx1 mRNA at this time [34].
Results support the present concept that the peak signal of IFNτ from the developing conceptus occurred earlier in MIX supplemented heifers, as there was a stark decrease in the abundance of mRNA transcripts encoding IRF9.Notably, Irf9 mRNA is typically more abundant between days 14-18 compared to days 25-40 of pregnancy in cattle, coincident with decreases in mRNA encoding IRF3, MX1α, MX1β, and MX2 in CAR at these two time points [59].To corroborate these findings, quantifying IFNτ in the uterine fluid would provide an indication of the concentration of this protein in the histotroph.Unfortunately, in the present study, recovery of the uterine fluid at sample collection proved too inconsistent for subsequent quantification of the concentrations of IFNτ protein.
For successful attachment and implantation to occur, there must also be restructuring of the endometrium and evasion of the maternal immune system by the developing conceptus [36,37].Functional analysis of the significantly affected mRNA transcripts in the present study identified significant changes in the "spliceosomal cycle", "cell death and survival", "cellular assembly and organization", and "autophagy" pathways, suggesting that significant structural and functional changes are occurring.The canonical pathway most predicted to be upregulated in MIX was the spliceosomal cycle.Most genes are transcribed as pre-mRNAs containing introns that must be spliced out to form the mRNA for translation [60][61][62].Splicing of pre-mRNA is facilitated by the enzymatic activity of the spliceosome, consisting of five ribonucleoprotein subunits (snRNPs) and various protein co-factors [63,64].This cycle is dependent on ATP [63] and is heavily regulated to ensure accurate splicing [60].Although the spliceosome cycle is significantly affected by the form of supplemental Se, the physiological relevance of this effect is unclear at this time.
Effective attachment and implantation also require a transient evasion of the maternal immune system to allow the fetal allograft to closely interact with the maternal interface [38].The control of this mechanism appears to be modulated predominantly by the endometrium, although it may also be stimulated by IFNτ and other signals from the developing conceptus [37,65].This immunotolerant environment is the result of a careful balance between compositional changes in the innate immune system consisting of granulocytes, monocytes, and dendritic cells, which may remain constitutively responsive to protect against pathogens [38], and the adaptive immune system, which must be repressed so as to not reject the invading conceptus [37,38].Ligands from the conceptus as well as prostaglandins act on the endometrium that has been primed by P4 to effect changes in both the innate and adaptive immune responses.There is a significant increase in the abundance of the MHC II-expressing molecules, macrophages, and dendritic cells, as well as NK cells that migrate into the uterine epithelial, endothelial, and stomal cells in the presence of elevated levels of P4.Additionally, there is an increased abundance of molecules that promote immune tolerance such as indoleamine 2,3-dioxygenase 1 (IDO1) and, subsequently, kynurenine, IL10, TGFB, and various cell surface receptors (CTLA4, PD-L1, and LAG3).Furthermore, there is increased expression of several ISGs associated with immune tolerance in the uterus and circulating immune cells [37].
Of interest, IDO1 catalyzes the rate-limiting step in the catabolism of tryptophan via the kynurenine pathway, although it is also critical in immune tolerance and protecting the embryo from rejection [37,66].The relative abundance of this transcript was significantly reduced in MIX versus ISe CAR samples.When comparing pregnancy to cyclic cows, IDO1 is significantly more abundant on day 17 of gestation [66,67], but then dramatically declines by day 20 [67].This elevated expression of IDO1 is coupled with a decrease in the relative concentration of tryptophan: kynurenine [66], and kynurenine can further affect immune tolerance [68].Relatedly, mRNA encoding ISG20 and SCARA1 were significantly lower, and mRNA encoding the transforming growth factor beta 1 (TGFB1) tended to be lower in the MIX-supplemented heifers, which is indicative of a decreased antiviral response and an increased tolerance for cell invasion [37,38].
In combination, these findings suggest more advanced preparation of the developing conceptus to evade immune rejection by the maternal endometrium.This supports our hypothesis that there is a shift in the timing of peak IFNτ signaling (MRP), and advancement of the conceptus and endometrium to better support impending attachment, implantation, and placentation.

Materials and Methods
This project (protocol number 2017-2828) was approved by the University of Kentucky's Institutional Animal Care and Use Committee.

Animals and Experimental Procedure
Angus-cross heifers (N = 20) received 45 days of Se depletion with no supplemental Se, followed by 45 days of Se repletion with a vitamin-mineral mix containing 35-ppm Se as only ISe to reestablish systemic Se to an adequate status [69,70].The heifers were then randomly assigned to supplemental Se treatment formulated with a vitamin-mineral mix containing 35 ppm Se as either inorganic Se (n = 10, ISe, sodium selenite, Prince Agri Products, Inc., Quincy, IL, USA) or a 1:1 ratio of ISe and OSe (n = 10, MIX, SEL-PLEX; Alltech, Inc., Nicholasville, KY, USA).Heifers received dietary treatments for at least 90 days prior to estrous synchronization.
Whole blood was collected via jugular venipuncture prior to depletion, repletion, and the start of treatment.It was also recovered bimonthly during treatment until the end of experimentation to confirm Se adequate status throughout.Whole blood was also collected from each heifer at estrus and prior to slaughter (d 17).Total blood Se was quantified by the University of Kentucky's Veterinary Diagnostics Laboratory (Lexington, KY, USA) using an Agilent 7900 inductively coupled plasma-mass spectrometer [71].Throughout experimentation, all animals maintained Se-adequate status [69,70], and there tended to be an effect of Se treatment (p = 0.07) as previously described in Crites et al. [34].

Experimental Regimen and Serum Collection
After at least 90 days of Se supplementation with each respective treatment (ISe vs. MIX), heifers were randomly injected with one or two doses of dinoprost tromethamine (25 mg, Lutalyse, Zoetis, Parsippany, NJ, USA) to induce regression of the CL.Heifers were then monitored daily for behavioral estrus (d 0) by visual observation as well as using CowManager technology (Gerverscop 9, Harmelen, The Netherlands) to predict the timing of estrus.The presence of a preovulatory follicle (estrus, d 0) was confirmed via transrectal ultrasonography using a 5-8 MHz linear transducer (Ibex Pro, E.I. Medical Imaging, Loveland, CO, USA) prior to artificial insemination at 0 h, 12 h, and 24 h using commercially available frozen semen from a single established high-fertility bull.
The reproductive tract (ovaries, uterus, and cervix) was collected from each heifer on d 17 following euthanasia via captive bolt stunning and exsanguination at the USDA inspected University of Kentucky Meat Laboratory.Initially, the uterus was flushed to collect the conceptus as described in Crites et al. [34].An intact conceptus was recovered from six heifers in each treatment group (ISe, n = 6; MIX, n = 6), and only these heifers were used for analyses in the present study.
Subsequently, caruncular endometrial samples were collected from the uterine horn ipsilateral to the ovary bearing the CL.This respective uterine horn was cut longitudinally to expose the uterine lumen, and an 8 mm biopsy punch (Integra LifeSciences Production Corporation, Mansfield, MA, USA) was used to collect CAR endometrial samples from all animals.These samples were flash-frozen in liquid N 2 and stored at −80 • C to be used for RNA extraction, RNA sequencing, and for real-time PCR (qPCR) analysis.

RNA Extraction
Total RNA was extracted from approximately 200 mg CAR samples of the endometrium using TRIzol Reagent (Invitrogen Corporation, Carlsbad, CA, USA).All samples had high purity, with 260/280 absorbance ratios of 1.88 or greater as quantified using a NanoDrop ND-100 Spectrophotometer (NanoDrop Technologies, Wilmington, DE, USA).

RNA-Sequencing
Transcriptomic analysis using RNA-sequencing was conducted using mRNA extracted from CAR tissue.Library preparation was performed by Zymo Research Corporation (Irvine, CA, USA).Initially, 500 ng of total RNA was used to construct the total RNA-Seq libraries.The method described in Bogdanova et al. [72] was used to remove ribosomal RNA (rRNA) with some modifications, and libraries were prepared using the Zymo-Seq RiboFree Total RNA Library Prep Kit (Zymo Research Corporation, Irvine, CA, USA).The RNA-Seq libraries were then sequenced using an Illumina NovaSeq with a sequencing depth of at least 30 million read pairs per sample.
Count data were uploaded to Integrated Differential Expression and Pathway Analysis (iDEP.96,[82]), and lowly expressed transcripts were removed at <1.0 count per million (CPM).Data were log2 transformed and then subjected to principal component analysis (PCA) and hierarchical clustering of the top 1000 expressed genes to visualize sample variation.For all samples, the average mapping percentage of the identified transcripts was 92.33%.The total counts for all samples are displayed in Figure 4A, and data for each sample were normally distributed (Figure 4B).The average correlation among all samples was 0.97, as indicated in Figure 5. Differentially expressed gene/transcript (DEG) expression analysis was completed using DESeq2 v1.28.0 [83], which uses the Wald test for hypothesis testing.At the significance level p < 0.05, a total of 2029 gene transcripts demonstrated an effect of treatment with a false discovery rate (FDR) of <40%.Presently, the high FDR appears to be due to the small sample size, sample variation, and relatively low proportion of DEGs in the MIX compared to ISe treatment groups.However, it is generally accepted that corroborating multiple gene transcripts after high-throughput analysis by qPCR provides validation of the observed changes [84].
Count data were uploaded to Integrated Differential Expression and Pathway Analysis (iDEP.96,[82]), and lowly expressed transcripts were removed at <1.0 count per million (CPM).Data were log2 transformed and then subjected to principal component analysis (PCA) and hierarchical clustering of the top 1000 expressed genes to visualize sample variation.For all samples, the average mapping percentage of the identified transcripts was 92.33%.The total counts for all samples are displayed in Figure 4A, and data for each sample were normally distributed (Figure 4B).The average correlation among all samples was 0.97, as indicated in Figure 5. Differentially expressed gene/transcript (DEG) expression analysis was completed using DESeq2 v1.28.0 [83], which uses the Wald test for hypothesis testing.At the significance level p < 0.05, a total of 2029 gene transcripts demonstrated an effect of treatment with a false discovery rate (FDR) of <40%.Presently, the high FDR appears to be due to the small sample size, sample variation, and relatively low proportion of DEGs in the MIX compared to ISe treatment groups.However, it is generally accepted that corroborating multiple gene transcripts after high-throughput analysis by qPCR provides validation of the observed changes [84].To assess global effects of Se treatment on the abundance of CAR transcripts at MRP, DEGs identified in DESeq2 analysis were analyzed for canonical, functional, and network analyses using QIAGEN's Ingenuity Pathway Analysis (IPA, QIAGEN, Redwood City, CA, USA, http://www.ingenuity.com(accessed on 10/03/2023)).The raw data (FASTQ files) of this manuscript were deposited into the National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA, accession number PRJNA1029628, release date: 30 April 2024).

Real-Time PCR Analysis
Real-time PCR (qPCR) was used to quantify the relative abundance of mRNA in CAR samples for identified genes that differed according to the RNA-Seq analysis and to characterize mRNA transcript abundance of selenoproteins and related receptors in CAR at MRP using a technique that is routinely used in our laboratory [10,34,43].Approximately 1 µg of RNA from each sample was reverse-transcribed into cDNA using SuperScript IV VILO Master Mix with ezDNAse Enzyme (Invitrogen by Thermo Fisher Scientific, Vilnius, Lithuania).Each sample was compared to a no-reverse-transcription control to ensure that subsequent qPCR results were not a result of genomic DNA contamination.
Subsequently, the relative abundance of mRNA transcripts Ifit2, Irf9, Isg20, Scara5, Scarb1, Timp2, Timp3, and Trim56 was determined to corroborate RNA-Seq results herein.Additionally, corroborating qPCR for the mRNA transcripts Ifit3, Irf1, Isg15, Oas2, and Rsad2 were published previously in Crites et al. [34] in a targeted mRNA analysis prior to conducting RNA sequencing in the present study.To assess global effects of Se treatment on the abundance of CAR transcripts at MRP, DEGs identified in DESeq2 analysis were analyzed for canonical, functional, and network analyses using QIAGEN's Ingenuity Pathway Analysis (IPA, QIAGEN, Redwood City, CA, USA, http://www.ingenuity.com(accessed on 10/03/2023)).The raw data (FASTQ files) of this manuscript were deposited into the National Center for Biotechnology Information (NCBI) Sequence Read Archive (SRA, accession number PRJNA1029628, release date: 30 April 2024).

Real-Time PCR Analysis
Real-time PCR (qPCR) was used to quantify the relative abundance of mRNA in CAR samples for identified genes that differed according to the RNA-Seq analysis and to characterize mRNA transcript abundance of selenoproteins and related receptors in CAR at MRP using a technique that is routinely used in our laboratory [10,34,43].Approximately Primers were designed against each respective RefSeq using the NCBI Primer-BLAST tool (https://www.ncbi.nlm.nih.gov/tools/primer-blast/(accessed on 9 June 2023)).Target products in cDNA were verified by DNA sequencing at ACGT, Inc. (Wheeling, IL, USA), and sequencing results were compared to each respective primer template using the NCBI Nucleotide-BLAST tool (https://blast.ncbi.nlm.nih.gov/Blast.cgi?PROGRAM= blastn&BLAST_SPEC=GeoBlast&PAGE_TYPE=BlastSearch (accessed on 8 August 2023)).The GenBank accession numbers, forward and reverse primer sequences, amplicon length of each product, and product identify for each transcript of interest are listed in Appendix A. To perform qPCR analysis, a total volume of 25 µL containing 5 µL of cDNA, 1 µL of a 10 µM stock of each primer (forward and reverse), 12.5 µL of 2 × SYBR Green PCR Master Mix (iTaq Universal SYBR Green Supermix, BIO-RAD, Hercules, CA, USA), and 5.5 µL of nuclease-free water was used.Reactions were conducted using a Bio-Rad CFX Maestro thermal cycler (Bio-Rad, Hercules, CA, USA).
The relative abundance of each transcript was quantified using the 2 −∆∆CT method [85] using transcripts for three consecutively expressed normally distributed genes that were not affected by Se treatment to normalize the data.For selenoprotein and selenoprotein P mRNA analyses, Gapdh, Hprt1, and Sdha transcripts were used as reference genes.Alternatively, for transcriptomics corroboration mRNA analysis, Gapdh, Hprt1, and RPS11 transcripts were the respective reference genes.Data were normalized to the expression level of ISe.A total of six heifers per treatment were analyzed via qPCR, and all reactions were performed in triplicate.

Statistical Analysis
Data are presented as least square means (±SEM), with individual heifer as the experimental unit.Data were analyzed for normal distribution and homogeneity.When appropriate, qPCR data were natural log transformed given that the samples were not distributed normally.For selenoprotein mRNA analysis, the transcripts Dio2, Gpx2, Selh, Selr, Selw, and Tfrc were natural log transformed, and for mRNA analysis to corroborate transcriptomic data, Isg20 was natural log transformed due to being non-normally distributed.The relative expression for each transcript identified in RNA-sequencing results was subjected to the Wald test using DESeq2, as described in Section 4.4."RNA-Sequencing" above.
To determine the effect of form of Se on concentrations of each mRNA transcript, data were analyzed using Student's t-test with the PROC TTest procedure of SAS statistical software package (version 9.4; SAS Institute, Inc., Cary, NC, USA), and each treatment group contained six heifers (ISe, n = 6; MIX, n = 6).Results were considered statistically significant at p ≤ 0.05 or with a tendency to differ at 0.05 < p ≤ 0.10.

Conclusions
This research had two distinct objectives.Initially, we performed a targeted mRNA analysis to quantify the relative abundance of mRNA encoding selenoproteins and selenoprotein P receptors in the CAR at d 17 of gestation (MRP) in response to ISe-or MIXsupplemental Se.Subsequently, we employed a global transcriptomic approach to expand upon previous findings in our lab that identified key changes in the abundance of mRNA encoding the P4 and IFNτ-induced DGAT2, FGF2, IFIT3, ISG15, MX1, OAS2, and RSAD2 in response to MIX-form Se diet compared to ISe, which occurred synonymously with longer conceptuses in MIX-supplemented heifers [34].Therefore, the second objective of this study was to examine transcriptomic changes using RNA sequencing analysis in the CAR tissue of the endometrium at MRP, testing the hypothesis that there will be changes in the endometrium indicating an advancement in CAR preparation for implantation of the conceptus.Results from the current transcriptomic analysis correspond with this hypothesis, as evidenced by the significant downregulation of interferon signaling through the JAK/STAT1/2 signal transduction pathway and downregulation of mRNA encoding the classical ISGs: IFIT1, IFIT2, IFIT3, IRF1, IRF9, ISG15, ISG20, OAS2, and RSAD2.This presumed shift in the timing of MRP appears plausible given the additional findings of changes in the cellular organization and function, and immune evasion occurring in MIX-supplemented heifers compared to those supplemented with ISe alone. 1 These contents are associated with each gene symbol and are the accession numbers of the sequences retrieved from the NCBI RefSeq database for designing primers and probes. 2All qPCR products were validated by sequencing.The identity values (%) presented are the base pair ratios between the total amplicon length and the number of identical base pairs. 1 These contents are associated with each gene symbol and are the accession numbers of the sequences retrieved from the NCBI RefSeq database for designing primers and probes. 2All qPCR products were validated by sequencing.The identity values (%) presented are the base pair ratios between the total amplicon length and the number of identical base pairs.* Primers are reported in [34].

Figure 1 .
Figure 1.(A) Score plot and (B) hierarchical matrix of data derived from RNA-Seq analys sample of caruncular (CAR) tissue from heifers supplemented with a vitamin-mineral mi ing 35 ppm as ISe (n = 6) or 1:1 mixture of ISe:OSe (MIX, n = 6).(A) Principal componen accounted for 25% of the variance, and principal component 2 (PC#2) accounted for 12% o iation among samples.(B) Hierarchical cluster analysis was composed of the top 1000 diff

Figure 1 .
Figure 1.(A) Score plot and (B) hierarchical matrix of data derived from RNA-Seq analysis of each sample of caruncular (CAR) tissue from heifers supplemented with a vitamin-mineral mix containing 35 ppm as ISe (n = 6) or 1:1 mixture of ISe:OSe (MIX, n = 6).(A) Principal component 1 (PC#1) accounted for 25% of the variance, and principal component 2 (PC#2) accounted for 12% of the variation among samples.(B) Hierarchical cluster analysis was composed of the top 1000 differentially expressed transcripts using iDEP.96.Significance was determined by Wald test, with significance declared at p < 0.05 for all.

Figure 2 .
Figure 2. Top IPA-identified canonical pathways of genes differentially expressed from CAR of heifers supplemented a vitamin-mineral mix containing 35 ppm as ISe (n = 6) or 1:1 mixture of ISe:OSe (MIX, n = 6).Pathways are ranked by z-score, with orange indicating a pathway that is predicted to be upregulated (positive z-score) and the blue indicating a pathway that is predicted to be downregulated (negative z-score).

Figure 2 .
Figure 2. Top IPA-identified canonical pathways of genes differentially expressed from CAR of heifers supplemented a vitamin-mineral mix containing 35 ppm as ISe (n = 6) or 1:1 mixture of ISe:OSe (MIX, n = 6).Pathways are ranked by z-score, with orange indicating a pathway that is predicted to be upregulated (positive z-score) and the blue indicating a pathway that is predicted to be downregulated (negative z-score).

Figure 3 .
Figure 3. Mechanistic illustration of downregulation of interferon-responsive genes an JAK/STAT1/2 signaling pathway.Red arrows demonstrate lower transcript abundance in MIX-Se compared to ISe-supplemented heifers.All transcripts included were significant via transcriptomi analysis, and an asterisk (*) indicates those transcripts that were corroborated via qPCR at p ≤ 0.05 Results of qPCR analysis for MX1 are reported in Crites et al. [34].

Figure 3 .
Figure 3.Mechanistic illustration of downregulation of interferon-responsive genes and JAK/STAT1/2 signaling pathway.Red arrows demonstrate lower transcript abundance in MIX-Secompared to ISe-supplemented heifers.All transcripts included were significant via transcriptomic analysis, and an asterisk (*) indicates those transcripts that were corroborated via qPCR at p ≤ 0.05.Results of qPCR analysis for MX1 are reported in Crites et al.[34].

Figure 4 .
Figure 4.The (A) total read counts and (B) distribution of log2 transformed data derived from RNA-Seq analysis of each sample of caruncular (CAR) tissue from heifers supplemented with a vitaminmineral mix containing 35 ppm as ISe (n = 6) or 1:1 mixture of ISe:OSe (MIX, n = 6).

Figure 4 .
Figure 4.The (A) total read counts and (B) distribution of log2 transformed data derived from RNA-Seq analysis of each sample of caruncular (CAR) tissue from heifers supplemented with a vitamin-mineral mix containing 35 ppm as ISe (n = 6) or 1:1 mixture of ISe:OSe (MIX, n = 6).

27 Figure 5 .
Figure 5. Correlation matrix of data derived from RNA-Seq analysis of each sample of caruncular (CAR) tissue from heifers supplemented with a vitamin-mineral mix containing 35 ppm as ISe (n = 6) or 1:1 mixture of ISe:OSe (MIX, n = 6).The average correlation among all samples was 0.97.

Figure 5 .
Figure 5. Correlation matrix of data derived from RNA-Seq analysis of each sample of caruncular (CAR) tissue from heifers supplemented with a vitamin-mineral mix containing 35 ppm as ISe (n = 6) or 1:1 mixture of ISe:OSe (MIX, n = 6).The average correlation among all samples was 0.97.
Int. J. Mol.Sci.2023, 24, x FOR PEER REVIEW MIX) of ISe and OSe (1:1 sodium selenite and SEL-PLEX).Selenium was supplemente ment groups ad libitum.2DataforDio2,Gpx2, Selh, Selr, Selw, and Tfrc were natural log tra due to having not-normal distribution.3Dataare expressed as a ratio of MIX relative t values are associated with Student's t-test using the PROC TTEST procedure of SAS statis ware package (version 9.4, SAS Institute Inc., Cary, NC, USA) at n = 6 per treatment.
2For statistical analysis, p-values were determined using the Wald test in DESeq2.

Table 1 .
Primer sets and product identities of qPCR analysis of selenoprotein and selenoprotein P transcripts.

Table 2 .
Primer sets and product identities of qPCR analysis of IFNτ-responsive transcripts and all reference transcripts.