Advertisement

Blood transcriptome analysis and identification of genes associated with supernumerary teats in Chinese Holstein cows

  • Author Footnotes
    * These authors contributed equally to this work.
    Q.Z. Chen
    Footnotes
    * These authors contributed equally to this work.
    Affiliations
    Laboratory of Animal Genetics, Breeding and Reproduction, Ministry of Agriculture of China, National Engineering Laboratory of Animal Breeding, College of Animal Science and Technology, China Agricultural University, 100193, Beijing, China
    Search for articles by this author
  • Author Footnotes
    * These authors contributed equally to this work.
    M.Y. Yang
    Footnotes
    * These authors contributed equally to this work.
    Affiliations
    Laboratory of Animal Genetics, Breeding and Reproduction, Ministry of Agriculture of China, National Engineering Laboratory of Animal Breeding, College of Animal Science and Technology, China Agricultural University, 100193, Beijing, China
    Search for articles by this author
  • X.Q. Liu
    Affiliations
    Laboratory of Animal Genetics, Breeding and Reproduction, Ministry of Agriculture of China, National Engineering Laboratory of Animal Breeding, College of Animal Science and Technology, China Agricultural University, 100193, Beijing, China
    Search for articles by this author
  • J.N. Zhang
    Affiliations
    Laboratory of Animal Genetics, Breeding and Reproduction, Ministry of Agriculture of China, National Engineering Laboratory of Animal Breeding, College of Animal Science and Technology, China Agricultural University, 100193, Beijing, China
    Search for articles by this author
  • S.Y. Mi
    Affiliations
    Laboratory of Animal Genetics, Breeding and Reproduction, Ministry of Agriculture of China, National Engineering Laboratory of Animal Breeding, College of Animal Science and Technology, China Agricultural University, 100193, Beijing, China
    Search for articles by this author
  • Y.J. Wang
    Affiliations
    Laboratory of Animal Genetics, Breeding and Reproduction, Ministry of Agriculture of China, National Engineering Laboratory of Animal Breeding, College of Animal Science and Technology, China Agricultural University, 100193, Beijing, China
    Search for articles by this author
  • W. Xiao
    Affiliations
    Beijing Animal Husbandry Station, Beijing, 100029, China
    Search for articles by this author
  • Y. Yu
    Correspondence
    Corresponding author
    Affiliations
    Laboratory of Animal Genetics, Breeding and Reproduction, Ministry of Agriculture of China, National Engineering Laboratory of Animal Breeding, College of Animal Science and Technology, China Agricultural University, 100193, Beijing, China
    Search for articles by this author
  • Author Footnotes
    * These authors contributed equally to this work.
Open AccessPublished:October 11, 2022DOI:https://doi.org/10.3168/jds.2022-22346

      ABSTRACT

      In dairy cows, supernumerary teats (SNT) are not desired as they are considered a repository for bacteria; thus, SNT are a risk factor for mastitis. Supernumerary teats are a heritable oligo- or polygenic trait. The incidence of SNT in offspring must be reduced by genomic selection. However, in modern dairy farming, farmers often ignore the effects of SNT on cows. The study aimed to elucidate the effects of SNT on dairy cows from the blood transcriptome level and identify genes associated with SNT in Chinese Holstein cows. We selected 6 SNT cows (Yes) and 6 non-SNT cows (No). In the 6 SNT cows, 3 cows had 1 SNT (One) and 3 cows had 2 SNT (Two). They were divided into 3 comparison groups (One vs. No; Two vs. No; and Yes vs. No). RNA was extracted from blood white membrane cells of 12 cows, and RNA sequencing was performed. Differential gene expression analysis based on the negative binomial distribution was used to detect differentially expressed genes in the One versus No and Two versus No comparison groups. Genes that were significantly upregulated or downregulated both in the One versus No and Two versus No groups (shared genes, SG) were obtained for further analysis. We also performed gene set enrichment analysis for all genes expressed in the Yes versus No group, correlation analysis between SG and the hematological parameters, protein–protein interaction network analysis of SG to select hub genes, and alternative splicing analysis for Yes versus No group to explore the functions of differentially spliced genes. We detected 289 SG. Gene set enrichment analysis, gene ontology, and the Kyoto Encyclopedia of Genes and Genomes enrichment analysis results showed that SNT affect immunity, inflammation, and lactation-related pathways in dairy cows. Correlation analysis showed that LOC104968484, SLC25A6, GADD45G, BAX, APAF1, ATM, XIAP, MDM4, BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1 are associated with white blood cell count and absolute value of lymphocytes in SNT cows only, so they might be genes associated with SNT in Chinese Holstein cows. We found 2 genes (BAX and MDM4) were also differentially spliced genes. However, the causal relationship between these genes and the SNT phenotype needs to be further studied. This study is the first to reveal the adverse effects of SNT on dairy cows at a transcriptional level, and the genes we found can be used as a reference for further searching for candidate genes for the SNT phenotype.

      Key words

      INTRODUCTION

      Supernumerary teats (SNT) refers to teats in addition to normally developing functional teats (
      • Pausch H.
      • Jung S.
      • Edel C.
      • Emmerling R.
      • Krogmeier D.
      • Götz K.-U.
      • Fries R.
      Genome-wide association study uncovers four QTL predisposing to supernumerary teats in cattle.
      ;
      • Martin P.
      • Palhiere I.
      • Tosser-Klopp G.
      • Rupp R.
      Heritability and genome-wide association mapping for supernumerary teats in French Alpine and Saanen dairy goats.
      ). In humans, SNT are associated with renal abnormalities and urinary diseases (pyuria, hematuria, and urolithiasis;
      • Méhes K.
      Association of supernumerary nipples with other anomalies.
      ;
      • Hersh J.
      Association of supernumerary nipples and renal anomalies—Reply.
      ;
      • Brown J.
      • Schwartz R.A.
      Supernumerary nipples and renal malformations: A family study.
      ). In pigs, nonfunctional teats, which fail to produce milk and nurture young, are equivalent to SNT. The mammary glands of pigs are underdeveloped and cannot store milk, which is a limiting factor in piglet growth and survival (
      • Hurley W.L.
      Review: Mammary gland development in swine: Embryo to early lactation.
      ).
      • Zhang L.
      • Peng F.
      • Yu F.
      • Wan L.
      • Zhou Z.Q.
      Expression of ESR1, PRLR, GHR, and IGF1R in mammary glands of Hu sheep with four teats.
      compared the mammary tissues of Hu sheep with or without SNT at different stages. During pregnancy and lactation, the mammary alveoli were denser in Hu sheep with SNT than in sheep without SNT (
      • Zhang L.
      • Peng F.
      • Yu F.
      • Wan L.
      • Zhou Z.Q.
      Expression of ESR1, PRLR, GHR, and IGF1R in mammary glands of Hu sheep with four teats.
      ). The number of mammary alveoli is proportional to the lactation ability of the mammary gland (
      • Brisken C.
      • Rajaram R.D.
      Alveolar and lactogenic differentiation.
      ). Thus, Hu sheep with SNT may have higher milk yield than those without.
      However, cows with SNT do not fare well. First, SNT not only have the potential to interfere with milking (
      • Mackenzie K.J.J.
      • Marshall F.H.A.
      On the presence of supernumerary mammary glands in cows and on their functional activity.
      ;
      • Skjervold H.
      Supernumerary teats in cattle.
      ), but they are also considered a reservoir for bacteria; thus, SNT are a risk factor for mastitis (
      • Joerg H.
      • Meili C.
      • Ruprecht O.
      • Bangerter E.
      • Burren A.
      • Bigler A.
      A genome-wide association study reveals a QTL influencing caudal supernumerary teats in Holstein cattle.
      ;
      • Butty A.M.
      • Frischknecht M.
      • Gredler B.
      • Neuenschwander S.
      • Moll J.
      • Bieber A.
      • Baes C.F.
      • Seefried F.R.
      Genetic and genomic analysis of hyperthelia in Brown Swiss cattle.
      ). In addition, SNT affect the milk performance of cows.
      • Wen H.
      • Luo H.
      • Mi S.
      • Liu X.
      • Wang Y.
      • Xiao W.
      • Yu Y.
      Effect of supernumerary teat on milk performance and genetic expression analysis of LGR5 gene in dairy cows.
      compared the milk performance of 325 SNT cows and 737 non-SNT cows; specifically, their results showed that the daily milk yield and peak milk of SNT cows are lower than those of non-SNT cows, and the gap between them can reach 2 kg in summer (
      • Wen H.
      • Luo H.
      • Mi S.
      • Liu X.
      • Wang Y.
      • Xiao W.
      • Yu Y.
      Effect of supernumerary teat on milk performance and genetic expression analysis of LGR5 gene in dairy cows.
      ).
      • Yang M.
      • Wen H.
      • Mi S.
      • Yu Y.
      The influence of supernumerary teats on milk production performance and blood routine in Chinese Holsteins.
      also showed that the daily milk yield, peak milk, and pre-milk yield of SNT cows are significantly lower than those of non-SNT cows (
      • Yang M.
      • Wen H.
      • Mi S.
      • Yu Y.
      The influence of supernumerary teats on milk production performance and blood routine in Chinese Holsteins.
      ). The existence of SNT affects the economic benefit of dairy farming. In modern dairy farming, SNT are mainly treated by artificial cutting. However, this method can increase the stress response of cattle and easily cause wound infections due to improper handling. Artificial cutting fails to meet the requirements of modern animal health and welfare (
      • Fraser D.
      Toward a global perspective on farm animal welfare.
      ).
      • Palacios C.
      • Abecia J.-A.
      Supernumerary teat removal can be avoided in dairy sheep.
      found that the milk yield of lactating sheep, whose SNT were cut off, significantly decreased compared with that of sheep with SNT (
      • Palacios C.
      • Abecia J.-A.
      Supernumerary teat removal can be avoided in dairy sheep.
      ). In cattle, the inheritance of SNT is oligo- or polygenic, with a heritability between 15 and 60% (
      • Brka M.
      • Reinsch N.
      • Kalm E.
      Frequency and heritability of supernumerary teats in German Simmental and German Brown Swiss cows.
      ;
      • Pausch H.
      • Jung S.
      • Edel C.
      • Emmerling R.
      • Krogmeier D.
      • Götz K.-U.
      • Fries R.
      Genome-wide association study uncovers four QTL predisposing to supernumerary teats in cattle.
      ). The heritability of SNT is different among breeds. For German Holsteins, the prevalence of SNT was 15%, and the posterior mean for the heritability was 0.09 with standard error of 0.08 (
      • Brka M.
      • Reinsch N.
      • Junge W.
      • Kalm E.
      Frequency and heritability of supernumerary teats in German Holsteins (Germany).
      ). The incidence of SNT in Chinese Holstein cows was about 9.8%, and the heritability was about 0.22 in our previous study (
      • Wen H.
      • Luo H.
      • Yang M.
      • Augustino S.M.A.
      • Wang D.
      • Mi S.
      • Guo Y.
      • Zhang Y.
      • Xiao W.
      • Wang Y.
      • Yu Y.
      Genetic parameters and weighted single-step genome-wide association study for supernumerary teats in Holstein cattle.
      ), which was medium heritability. Therefore, the incidence of SNT in cows must be reduced through genomic selection.
      Studies on SNT mostly focused on Simmental and German Brown Swiss cows, and few studies worked on Holstein cows, which have stronger lactation performance than other cows (
      • Brka M.
      • Reinsch N.
      • Kalm E.
      Frequency and heritability of supernumerary teats in German Simmental and German Brown Swiss cows.
      ). Functional genes associated with SNT were mainly identified by GWAS analysis (
      • Pausch H.
      • Jung S.
      • Edel C.
      • Emmerling R.
      • Krogmeier D.
      • Götz K.-U.
      • Fries R.
      Genome-wide association study uncovers four QTL predisposing to supernumerary teats in cattle.
      ;
      • Joerg H.
      • Meili C.
      • Ruprecht O.
      • Bangerter E.
      • Burren A.
      • Bigler A.
      A genome-wide association study reveals a QTL influencing caudal supernumerary teats in Holstein cattle.
      ;
      • Butty A.M.
      • Frischknecht M.
      • Gredler B.
      • Neuenschwander S.
      • Moll J.
      • Bieber A.
      • Baes C.F.
      • Seefried F.R.
      Genetic and genomic analysis of hyperthelia in Brown Swiss cattle.
      ). Our previous study identified 10 SNT-associated genes in BTA4 by single-step GWAS analysis (
      • Wen H.
      • Luo H.
      • Yang M.
      • Augustino S.M.A.
      • Wang D.
      • Mi S.
      • Guo Y.
      • Zhang Y.
      • Xiao W.
      • Wang Y.
      • Yu Y.
      Genetic parameters and weighted single-step genome-wide association study for supernumerary teats in Holstein cattle.
      ). However, transcriptome sequencing (RNA-seq) has not been used to carry out relevant studies on cows with SNT. Blood is an important material in the process of disease research, and its composition changes represent the health status of animals. In addition, hematological parameters (HP) are easy to obtain. In view of these, blood RNA-seq was performed on Chinese Holstein dairy cows with and without SNT to explore the effects of SNT on cows and search for genes associated with SNT.

      MATERIALS AND METHODS

      Ethics Statement

      The experimental procedures were approved by the Animal Welfare Committee of China Agricultural University, Beijing, China. Animal welfare and handling procedures were strictly followed throughout the experiment.

      Collection of Animals and Blood Samples

      The experimental cows came from a dairy farm of Beijing Shounong Group (Beijing, China). A total of 12 Chinese Holstein cows were selected based on their parity, lactation stage, and calving age. Half of the 12 cows were non-SNT cows (No), selected from non-SNT cows. The rest were SNT cows [Yes; including 3 cows with one SNT (One) and 3 cows with 2 SNT (Two)], selected from SNT cows. The experiment was divided into 3 comparison groups (One vs. No, Two vs. No, and Yes vs. No). About 20 mL of tail venous blood was collected from each cow. Ten milliliters of tail venous blood was used for RNA-seq, and the remaining 10 mL of tail venous blood was sent to Jinhaikeyu Company for HP testing (Beijing, China).

      mRNA Extraction and Sequencing

      Blood was immediately centrifuged at 1,057 × g for 15 min at 4°C to isolate white membrane layer and stored in liquid nitrogen for RNA-seq. We extracted total RNA from blood white membrane layer (mainly consists of platelets, lymphocytes, monocytes, and granulocytes) of 12 dairy cows by using TRIzol reagent (Invitrogen). The concentration and integrity of RNA were determined by NANODROP 2000 spectrophotometer and Bioanalyzer 2100 system (Agilent Technologies), respectively. The concentration of the RNA samples ranged from 150 ng/μL to 800 ng/μL, and their RNA integrity numbers were greater than 6 (Supplemental Table S1; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ), which was sufficient for gene expression analysis (
      • Cánovas A.
      • Rincón G.
      • Bevilacqua C.
      • Islas-Trejo A.
      • Brenaut P.
      • Hovey R.C.
      • Boutinaud M.
      • Morgenthaler C.
      • VanKlompenberg M.K.
      • Martin P.
      • Medrano J.F.
      Comparison of five different RNA sources to examine the lactating bovine mammary gland transcriptome using RNA-sequencing.
      ;
      • Gallego Romero I.
      • Pai A.A.
      • Tung J.
      • Gilad Y.
      RNA-seq: Impact of RNA degradation on transcript quantification.
      ;
      • Davila J.I.
      • Fadra N.M.
      • Wang X.
      • McDonald A.M.
      • Nair A.A.
      • Crusan B.R.
      • Wu X.
      • Blommel J.H.
      • Jen J.
      • Rumilla K.M.
      • Jenkins R.B.
      • Aypar U.
      • Klee E.W.
      • Kipp B.R.
      • Halling K.C.
      Impact of RNA degradation on fusion detection by RNA-seq.
      ). Novogene Co. Ltd. (Beijing, China) performed the cDNA library preparation and RNA sequencing, generating paired-end reads at 150 bp length on the Illumina HiSeq 4000 platform.

      RNA-seq Data Processing and Analysis

      We conducted quality assessment for raw reads using FastQC (version 0.11.3; http://www.bioinformatics.babraham.ac.uk/projects/fastqc/). Subsequently, we filtered adapters and reads with low quality by using the NGSQCToolkit software (version 2.3; https://github.com/mjain-lab/NGSQCToolkit/releases/tag/v2.3) to obtain clean reads. The proportion of bases with low quality in reads >30% were considered reads with low quality. Bases with low quality were those with a mass less than 20. The SAM files were observed after Hisat2 software (version 2.1.0; https://daehwankimlab.github.io/hisat2/download/) was used to align with the bovine reference genome (http://ftp.ensembl.org/pub/release-97/gtf/bos_taurus/Bos_taurus.ARS-UCD1.2.97.gtf.gz). The details of sequencing data are shown in Supplemental Table S2 (https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). The SAM files were then converted to BAM files by SAMtools (version 1.16.1; https://github.com/samtools/samtools/releases/tag/1.16.1). Rsubread (version 2.10.5; http://bioconductor.org/packages/release/bioc/html/Rsubread.html) software was used to generate the read counts for each sample. Genes that were not expressed in all samples were deleted, and read counts were corrected to fragments per kilobase of transcript per million mapped reads (Supplemental Table S3; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). The DESeq2 package was used to find significantly differentially expressed genes (DEG) in the One versus No and Two versus No groups [|log2-fold change| > 1, adjusted P-value (P-adjust) < 0.05]. Given that our study focused on stable SNT genetic markers of cows. The “marker” in our study implied significantly DEG shared by group One versus No and group Two versus No in the same direction of difference (shared genes, SG). In other words, SG were upregulated in both One versus No and Two versus No groups (up-SG) and significantly downregulated in both One versus No and Two versus No groups (down-SG). Therefore, further analysis of SG is meaningful. Gene set enrichment analysis (GSEA) was performed using all genes in the Yes versus No group by GSEABase package (version 1.36.0; http://www.bioconductor.org/packages/3.4/bioc/html/GSEABase.html).
      To further search for genes associated with SNT of Chinese Holstein cows, we performed correlation analysis between the expression level of SG and HP (the blood routine test data shown in Supplemental Table S4; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ) of cows using the lm function of R (version 4.2.1; https://cran.r-project.org/). To understand the relationship between proteins encoded by SG, protein–protein interaction (PPI) network analysis was performed via Cytoscape software (version 3.9.1; https://cytoscape.org/). CytoHubba within Cytoscape (version 3.9.1) was then used for hub gene analysis, and maximal clique centrality was used as the rank criteria. To carry out alternative splicing analysis, STAR software (version 2.7.9; https://github.com/alexdobin/STAR/releases/tag/2.7.9a) was used for genome assembly to obtain BAM files. The rMATS (version 4.1.2) software (http://rnaseq-mats.sourceforge.net/download.html) was used for alternative splicing analysis. Five types of splicing events were recognized by rMATS (version 4.1.2), namely, skipped exon, alternative 5′ splice site, alternative 3′ splice site, mutually exclusive exons, and retained intron. The rMATS default parameters were used for the differential splicing analysis of the Yes versus No group. Transcripts with P < 0.05 were considered differentially spliced transcripts. Intersections between differentially spliced genes (DSG) of different splicing events and SG were shown by UpSetR package (version 1.4.0; https://cran.rstudio.com/web/packages/UpSetR/). The results of DSG were visualized by rmats2sashimiplot software (version 2.0.4; https://github.com/Xinglab/rmats2sashimiplot/releases/tag/v2.0.4). Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analyses were performed for SG and DSG using the ClusterProfiler R package. Gene Ontology terms and KEGG pathways with P < 0.05 were obtained.

      Validation of Sequencing Results by Real-Time Quantitative PCR

      PrimeScript RT reagent Kit with gDNA Eraser (Takara Biotechnology Co. Ltd.) was used to produce cDNA from qualified RNA. The bovine gene GAPDH was used as the reference gene, and 12 DEG were randomly selected. Primers were designed through online website Primer-Blast, and the primers' information is shown in Table 1. The primer sequences were synthesized by Tsingke Biotechnology Co., Ltd. (Beijing, China). Three technical replicates were prepared for each cDNA sample, and the mean cycle threshold (Ct) value was used to calculate the genes' relative expression by the 2−ΔΔCt method. The significant differences between groups were examined by t-test.
      Table 1The primer information for real-time quantitative PCR (F = forward, R = reverse)
      GenePrimer sequenceProduct length (bp)
      ADAMDEC1F: AGCCAAGATGCCTACTGCAA237
      R: GCTTGGATTCATTGAGTTTGGA
      ADIRFF: GCGGCTACCCTTAACCACTC189
      R: TGGGTAGTCTTGGCAACCTG
      C1QCF: GCCCCTCCTGGTACTAAACCT140
      R: CTGGCAGTCCGTCATACCCA
      CCL3F: ACACTCCGTCTCGCAGCATC291
      R: GCTCCAGGTCGGTGATGTATT
      CLEC1BF: AGTGACGTCTGGAGATGGGA166
      R: CGCCTGCTTTCCTTTCACAC
      FSCN1F: CAAGGACGCCAGATGCTACT288
      R: ATGTTGTAGGCACCGTCTCG
      HECTD4F: CCTTGCAGGGAGAGACATTGA124
      R: ATGAGGGCAGTTGTGACCTG
      LOC104968634F: CCTACTCTGTGGTTCTTGTCAGAGG218
      R: TGATGTCCTCACAGACAACCC
      LOC618541F: GCTCTGACTGCAGGATCTCA222
      R: GACTGTGAACCCACCTACCG
      PLAC8AF: ACCCAGTTGTTTCACAGCCA193
      R: GTCCTCATGGCGACACTTGA
      PRSS57F: CATCACGCACCCAGACTACC190
      R: CAGTTCCTCGAAGTCGGACA
      RETNF: GCAGCCTTAAAAAGGGAGCG222
      R: TCCTTACTGCCCCAGGAACT
      GAPDHF: GGTGCTGAGTATGTGGTGGA119
      R: GGCATTGCTGACAATCTTGA

      RESULTS

      Differential Expression Profile of Genes Between Cows With and Without SNT

      A total of 397 and 1,559 DEG were detected in the One versus No and Two versus No groups, respectively (|log2-fold change| > 1, P-adjust < 0.05; Supplemental Tables S5, S6; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). As shown in the clustering heatmap (Figure 1A), the difference in the expression patterns of DEG in the 2 comparison groups was obvious. In addition, a total of 20,303 genes expressed in the Yes versus No group were used for principal component analysis and hierarchical clustering. The results showed that SNT could classify the groups as left side and right side, which illustrated SNT as one of the main influencing factors of the differences in samples (Figure 1B). In DEG, 134 and 406 genes were significantly upregulated in the One versus No and Two versus No groups, respectively, and most of the DEG were significantly downregulated in the One versus No (263) and Two versus No (1,153) groups (Figure 1C). Among them, we found 84 genes that were significantly upregulated in the One versus No and Two versus No groups (left panel in Figure 1D), and 205 genes that were significantly downregulated in the One versus No and Two versus No groups (right panel in Figure 1D).
      Figure thumbnail gr1
      Figure 1Differentially expressed genes (DEG) between cows with supernumerary teats (SNT) and non-SNT cows. (A) Heatmaps of DEG in the One versus No and Two versus No groups. (B) Principal component analysis (PCA) of 20,303 genes expressed in SNT and non-SNT cows. (C) Volcano plot of DEG in the One versus No and Two versus No groups. (D) Genes that were significantly upregulated both in the Two versus No and One versus No groups, and significantly downregulated both in the Two versus No and One versus No groups. Up = significantly upregulated genes; down = significantly downregulated genes; No = non-SNT cows; Yes = SNT cows; One = cows with 1 SNT; Two = cows with 2 SNT; up-SG = genes that were significantly upregulated both in the One versus No and Two versus No groups; down-SG = genes that were significantly downregulated both in the One versus No and Two versus No groups. FDR = false discovery rate.

      Genes Expressed Between SNT Cows and Non-SNT Cows, or SG Enriched in Immunity or Inflammation-Related Pathways

      First, we used all genes expressed in the Yes versus No group for GSEA. The results showed that the SNT phenotype was positively correlated with inflammation-related pathways (P < 0.05), such as systemic lupus erythematosus (SLE), neutrophil extracellular trap (NET) formation, necroptosis, p53 signaling pathway, NF-κB signaling pathway, IL-17 signaling pathway, Staphylococcus aureus infection, NOD-like receptor signaling pathway, and TNF signaling pathway (Figure 2A). Moreover, the SNT phenotype was negatively correlated with the prolactin signaling pathway and carbohydrate digestion and absorption (Figure 2A). Shared genes were selected for functional enrichment analysis. Gene Ontology and KEGG enrichment analyses were carried for SG. A total of 80 notable GO terms were enriched (P < 0.05; Supplemental Table S7; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). The top 40 terms of biological process were chosen based on the P-value for visualizations (Figure 2B), and these terms were mainly associated with immune regulation (e.g., T cell cytokine production, T cell-mediated immunity, and regulation of immune response), heme metabolic process, and post-transcriptional gene silencing. Posttranscriptional gene silencing may be associated with downregulation of most genes in SNT cows. A total of 22 pathways were enriched by SG (P < 0.05; Supplemental Table S8; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ), which were mainly associated with immunity, inflammation, and apoptosis (e.g., SLE, NET formation, apoptosis, and p53 signaling pathway and NF-κB signaling pathway; Figure 2C). In summary, these results indicated that SNT is associated with immunity and inflammation of dairy cows. Supernumerary teats may have adverse effects on the health of cows and affect lactation performance.
      Figure thumbnail gr2
      Figure 2Functional enrichment analysis of genes. (A) Gene set enrichment analysis was performed for 20,303 genes expressed in the Yes versus No group. The supernumerary teat (SNT) phenotype was significantly positively correlated with inflammation-related pathways and significantly negatively correlated with lactation-related pathways. (B) Biological process (BP) of gene ontology (GO) analysis of shared genes (SG; top40). The terms highlighted in red are those associated with immunity. (C) The left shows the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways enriched by SG, whereas the genes shown on the right are the genes in the corresponding pathway. These genes were chosen based on correlation analysis, which found that these genes were strongly correlated with white blood cell count and absolute value of lymphocytes only in SNT cows. The y-axis is corrected read counts of genes. No = non-SNT cows; Yes = SNT cows.

      Strong Association Between Genes and HP in SNT Cows Only

      To further search for genes associated with SNT, correlation analysis between SG and HP was conducted. We found that LOC104968484, SLC25A6, GADD45G, and BAX were significantly, or tended to be significantly, negatively correlated with white blood cell (WBC) count and absolute value of lymphocytes (W-SCC) only in SNT cows (green line in Figure 3A), and they were significantly upregulated in SNT cows (Figure 2C). By contrast, the expression levels of APAF1, ATM, XIAP, and MDM4, which were significantly downregulated in SNT cows (Figure 2C), were significantly, or tended to be significantly, positively correlated with WBC count and W-SCC in SNT cows (green line in Figure 3B). However, the upregulated genes (LOC104968484, SLC25A6, GADD45G, and BAX) and downregulated (APAF1, ATM, XIAP, and MDM4) genes did not have such a relationship with WBC count and W-SCC in non-SNT cows (red line in Figures 3A and 3B). More interestingly, any significant correlation was not observed between the 8 genes (APAF1, ATM, XIAP, MDM4, LOC104968484, SLC25A6, GADD45G, and BAX) and red blood cell count and platelet count in SNT or non-SNT (Figures 4A and 4B). The results suggested that the 8 genes (APAF1, ATM, XIAP, MDM4, LOC104968484, SLC25A6, GADD45G, and BAX) were highly associated with immunity of SNT cows, and they could be considered important genes associated with SNT (IGR-SNT).
      Figure thumbnail gr3
      Figure 3Correlation between the expression level of important genes associated with supernumerary teats (IGR-SNT) and (a) white blood cell (WBC) count and (b) absolute value of lymphocytes (W-SCC) in cows with supernumerary teats (SNT) and non-SNT cows, respectively. (A) Correlation between the expression level of upregulated IGR-SNT (LOC104968484, SLC25A6, GADD45G, and BAX) and WBC count and W-SCC. The upregulated IGR-SNT were significantly negatively correlated (P < 0.05) or tended to be negatively correlated with WBC count and W-SCC only in SNT cows. (B) Correlation between the expression level of downregulated IGR-SNT (APAF1, ATM, XIAP, and MDM4) and WBC count and W-SCC. The downregulated IGR-SNT were significantly positively correlated (P < 0.05) or tended to be positively correlated with WBC count and W-SCC only in SNT cows. The x-axis is the corrected read counts of IGR-SNT.
      Figure thumbnail gr4
      Figure 4Correlation between the expression level of important genes associated with supernumerary teats (IGR-SNT) and (m) red blood cell (RBC) count and (n) platelet (PLT) count in cows with supernumerary teats (SNT) and non-SNT cows, respectively. (A) Correlation between the expression level of upregulated IGR-SNT (LOC104968484, SLC25A6, GADD45G, and BAX) and RBC count and PLT count. The upregulated IGR-SNT had no correlation with RBC count and PLT count in SNT cows or non-SNT cows. (B) Correlation between the expression level of downregulated IGR-SNT (APAF1, ATM, XIAP, and MDM4) and RBC count and PLT count. The downregulated IGR-SNT had no correlation with RBC count and PLT count in SNT cows or non-SNT cows. The x-axis is the corrected read counts of IGR-SNT.
      We conducted PPI network analysis for SG, and the results showed that highly interacting proteins, including SMG1, BDP1, CEP350, MED13, AFF4, PHIP, TAOK1, ATM, MAP3K2, PIK3C2A, TAOK1, BRWD3, and ESR1, were encoded by genes that were significantly downregulated in SNT cows (Figure 5A). CytoHubba was then used for screening hub genes, and the rank criteria was selected using the maximal clique centrality method. In our study, the top 10 genes were considered hub genes, and a high degree of interaction was found between the hub genes (Figure 5B). We found 6 hub genes (BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) were more correlated with WBC count and W-SCC in SNT cows than in non-SNT cows (Supplemental Figure S1A; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). A significant correlation was not observed between 6 hub genes (BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) and red blood cell count and platelet count in SNT or non-SNT (Supplemental Figure S1B;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). Thus, the 6 hub genes (BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) were also likely to be associated with the immunity of SNT cows.
      Figure thumbnail gr5
      Figure 5Protein–protein interaction (PPI) network analysis of shared genes (SG). (A) PPI network of SG. Green nodes = significantly downregulated genes in cows with supernumerary teats (SNT); red nodes = significantly upregulated genes in SNT cows. Nodes with black circles were important genes associated with SNT (IGR-SNT); nodes with red circles were hub genes. (B) The top 10 hub genes ranked by the maximal clique centrality method, which was in CytoHubba software (version 2.7.9).
      In brief, the IGR-SNT and the 6 hub genes had a strong association with WBC count and W-SCC in SNT cows only. They could be considered genes associated with SNT. The differences between IGR-SNT and the 6 hub genes were that in SNT cows, the r values (correlation coefficients) of correlation between IGR-SNT, WBC count, and W-SCC were generally higher than that of the 6 hub genes (right panel in Supplemental Figure S1C;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). The P-value showed the opposite trend (left panel in Supplemental Figure S1C;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). In addition, based on the result of KEGG enrichment analysis for SG, the 6 hub genes (BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) were not enriched in any prominent pathways, but the IGR-SNT were all enriched (Figure 2C; Supplemental Table S8; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ).

      Hub Genes of SG were Highly Correlated with IGR-SNT

      Hub genes play an important role in the PPI network. In other words, the regulation of other genes is often influenced by the hub gene. Further analysis revealed that the expression levels of the 6 hub genes (BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) in SNT cows were almost significantly positively correlated with the downregulated IGR-SNT (APAF1, ATM, XIAP, and MDM4; Figure 6A). In SNT cows, although the P-value of correlation between the 6 hub genes and upregulated IGR-SNT (LOC104968484, SLC25A6, GADD45G, and BAX) was almost above 0.05, they all showed a negative correlation trend (Supplemental Figure S2; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). Thus, a strong regulatory relationship was found between the 6 hub genes and the IGR-SNT, and the relationship between the 6 hub genes and the IGR-SNT is shown in Figure 6B.
      Figure thumbnail gr6
      Figure 6The correlation between the expression level of the 6 hub genes and the important genes associated with supernumerary teats (IGR-SNT). (A) Correlation between the expression level of hub genes (BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) and downregulated IGR-SNT (APAF1, ATM, XIAP, and MDM4) in cows with supernumerary teats (SNT). The y-axis shows the corrected read counts of downregulated IGR-SNT. The x-axis shows corrected read counts of hub genes. (B) Schematic of the relationship between the 6 hub genes (BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) and IGR-SNT (LOC104968484, SLC25A6, GADD45G, BAX, APAF1, ATM, XIAP, and MDM4) in SNT cows. Plus sign = positive correlation; minus sign = negative correlation; arrow pointing down = downregulated genes in SNT cows; arrow pointing up = upregulated genes in SNT cows; WBC = white blood cell; W-SCC = absolute value of lymphocytes.

      Alternative Splicing Analysis of Genes Between SNT Cows and Non-SNT Cows

      To explore the mechanism of alternative splicing in the formation of SNT, rMATS were used to examine DSG of different splicing events. A total of 3,060 DSG were obtained in the Yes versus No group (P < 0.05). More than half of the DSG splicing type were skipped exon. The results of differentially spliced transcripts of different splicing events are shown in Supplemental Tables S9 through S13 (https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). A total of 57 SG were DSG (Figure 7A). A KEGG enrichment analysis was carried out for DSG. A total of 79 pathways were enriched (P-adjust < 0.05; Supplemental Table S14; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ). The top 30 pathways were mainly associated with immunity and inflammation (Figure 7B). More importantly, 2 IGR-SNT (BAX and MDM4) were identified as DSG (Figure 7C and D).
      Figure thumbnail gr7
      Figure 7Alternative splicing analysis of genes between cows with supernumerary teats (SNT) and non-SNT cows. (A) Intersection between differentially spliced genes (DSG) of different splicing events and shared genes (SG). (B) Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment analysis of DSG (top 30). The pathways highlighted in red were associated with immunity or inflammation. (C) Alternative splicing visualization of the 2 genes (BAX and MDM4), which were both differentially expressed and differentially spliced genes. A3SS = alternative 3′ splice site; A5SS = alternative 5′ splice site; RI = retained intron; MXE = mutually exclusive exons; SE = skipped exons.

      Validation of RNA-seq Results by Real-Time Quantitative PCR

      To verify the results of DEG obtained by RNA-seq, we randomly selected 12 genes from the One versus No and Two versus No groups to perform real-time quantitative PCR. The comparison results of RNA-seq and real-time quantitative PCR showed the high reliability of RNA-seq (Supplemental Figure S3; https://doi.org/10.6084/m9.figshare.20186543.v6;
      • Chen Q.
      Supplemental tables and figures. figshare. Figure.
      ).

      DISCUSSION

      This work is the first study to reveal the adverse effects of SNT on dairy cows at transcriptional level and explore the genes associated with SNT in Chinese Holstein cows by using RNA-seq technology.
      Inflammation-related pathways were activated, and prolactin signaling pathway and carbohydrate digestion and absorption were inhibited in SNT cows. The GSEA results showed that the SNT phenotype was positively correlated with inflammation-related pathways (e.g., IL-17 signaling pathway, Staphylococcus aureus infection, TNF signaling pathway, and NF-κB signaling pathway), and negatively correlated with prolactin signaling pathway and carbohydrate digestion and absorption. Staphylococcus aureus is an important gram-positive bacterial pathogen that causes IMI, and is associated with various forms of subclinical and clinical mastitis in dairy cows (
      • Liu B.
      • Li Q.
      • Gong Z.
      • Zhao J.
      • Gu B.
      • Feng S.
      Staphylococcus aureus lipoproteins play crucial roles in inducing inflammatory responses and bacterial internalization into bovine mammary epithelial cells.
      ). Supernumerary teats may constitute a minor, but notable, risk factor for mastitis (
      • Vasileiou N.G.C.
      • Mavrogianni V.S.
      • Petinaki E.
      • Fthenakis G.C.
      Predisposing factors for bacterial mastitis in ewes.
      ), which was consistent with our results, which found that the S. aureus infection pathway was activated in SNT cows. Studies have shown that SNT cows have poor milk performance (
      • Wen H.
      • Luo H.
      • Mi S.
      • Liu X.
      • Wang Y.
      • Xiao W.
      • Yu Y.
      Effect of supernumerary teat on milk performance and genetic expression analysis of LGR5 gene in dairy cows.
      ;
      • Yang M.
      • Wen H.
      • Mi S.
      • Yu Y.
      The influence of supernumerary teats on milk production performance and blood routine in Chinese Holsteins.
      ). Our result found that the prolactin signaling pathway and carbohydrate digestion and absorption were inhibited in SNT cows. Prolactin plays an important role in lactation, which can induce breast acinar cells to produce milk protein (
      • Zhang Z.
      • Wei Q.
      • Zeng Y.
      • Jia X.
      • Su H.
      • Lin W.
      • Xing N.
      • Bai H.
      • He Y.
      • Wang Q.
      Effect of Hordei Fructus Germinatus on differential gene expression in the prolactin signaling pathway in the mammary gland of lactating rats.
      ). Long-term inhibition of prolactin reduces milk yield of dairy cows (
      • Lacasse P.
      • Lollivier V.
      • Bruckmaier R.M.
      • Boisclair Y.R.
      • Wagner G.F.
      • Boutinaud M.
      Effect of the prolactin-release inhibitor quinagolide on lactating dairy cows.
      ;
      • Lollivier V.
      • Lacasse P.
      • Angulo Arizala J.
      • Lamberton P.
      • Wiart S.
      • Portanguen J.
      • Bruckmaier R.
      • Boutinaud M.
      In vivo inhibition followed by exogenous supplementation demonstrates galactopoietic effects of prolactin on mammary tissue and milk production in dairy cows.
      ). Carbohydrate is the main component of ruminant diet and can be divided into fibrous carbohydrate and NFC. Non-fibrous carbohydrate is an important dietary energy source for ruminants, which can affect the performance of ruminants (
      • Batajoo K.K.
      • Shaver R.D.
      Impact of nonfiber carbohydrate on intake, digestion, and milk-production by dairy-cows.
      ;
      • Gao X.
      • Oba M.
      Effect of increasing dietary nonfiber carbohydrate with starch, sucrose, or lactose on rumen fermentation and productivity of lactating dairy cows.
      ).
      • Batajoo K.K.
      • Shaver R.D.
      Impact of nonfiber carbohydrate on intake, digestion, and milk-production by dairy-cows.
      found that decreasing dietary NFC decreased DMI, milk protein percentage, and production in dairy cows. Thus, we speculated that the decreased lactation performance of SNT cows in other studies (
      • Wen H.
      • Luo H.
      • Mi S.
      • Liu X.
      • Wang Y.
      • Xiao W.
      • Yu Y.
      Effect of supernumerary teat on milk performance and genetic expression analysis of LGR5 gene in dairy cows.
      ;
      • Yang M.
      • Wen H.
      • Mi S.
      • Yu Y.
      The influence of supernumerary teats on milk production performance and blood routine in Chinese Holsteins.
      ) may be partly associated with the inhibition of prolactin and carbohydrate digestion and absorption in cows with SNT. Shared genes were used for functional enrichment analysis. Multiple immunity-related GO terms and inflammation-related pathways were enriched by SG, such as NET formation. Neutrophil extracellular traps are extracellular DNA fibers composed of histones and neutrophil antimicrobial proteins; particularly, they are released from neutrophils to kill pathogens (
      • Klopf J.
      • Brostjan C.
      • Eilenberg W.
      • Neumayer C.
      Neutrophil extracellular traps and their implications in cardiovascular and inflammatory disease.
      ). Although NET play an important role in host defense under infection, NET formation has recently been implicated in the pathophysiology of autoimmunity (
      • Grayson P.C.
      • Kaplan M.J.
      At the bench: Neutrophil extracellular traps (NETs) highlight novel aspects of innate immune system involvement in autoimmune diseases.
      ). For example, extensive NET formation has been shown to aggravate SLE (
      • Villanueva E.
      • Yalavarthi S.
      • Berthier C.C.
      • Hodgin J.B.
      • Khandpur R.
      • Lin A.M.
      • Rubin C.J.
      • Zhao W.
      • Olsen S.H.
      • Klinker M.
      • Shealy D.
      • Denny M.F.
      • Plumas J.
      • Chaperot L.
      • Kretzler M.
      • Bruce A.T.
      • Kaplan M.J.
      Netting neutrophils induce endothelial damage, infiltrate tissues, and expose immunostimulatory molecules in systemic lupus erythematosus.
      ). In our study, SLE was the most prominent pathway in KEGG enrichment. Furthermore, GSEA results revealed that NET formation and SLE were positively correlated with the SNT phenotype.
      • Pausch H.
      • Jung S.
      • Edel C.
      • Emmerling R.
      • Krogmeier D.
      • Götz K.-U.
      • Fries R.
      Genome-wide association study uncovers four QTL predisposing to supernumerary teats in cattle.
      and
      • Peng W.F.
      • Xu S.S.
      • Ren X.
      • Lv F.H.
      • Xie X.L.
      • Zhao Y.X.
      • Zhang M.
      • Shen Z.Q.
      • Ren Y.L.
      • Gao L.
      • Shen M.
      • Kantanen J.
      • Li M.H.
      A genome-wide association study reveals candidate genes for the supernumerary nipple phenotype in sheep (Ovis aries).
      both detected the Wnt signaling pathway, which plays an important role in the formation of SNT through GWAS analysis (
      • Pausch H.
      • Jung S.
      • Edel C.
      • Emmerling R.
      • Krogmeier D.
      • Götz K.-U.
      • Fries R.
      Genome-wide association study uncovers four QTL predisposing to supernumerary teats in cattle.
      ;
      • Peng W.F.
      • Xu S.S.
      • Ren X.
      • Lv F.H.
      • Xie X.L.
      • Zhao Y.X.
      • Zhang M.
      • Shen Z.Q.
      • Ren Y.L.
      • Gao L.
      • Shen M.
      • Kantanen J.
      • Li M.H.
      A genome-wide association study reveals candidate genes for the supernumerary nipple phenotype in sheep (Ovis aries).
      ). However, the Wnt pathway was not enriched by RNA-seq in our study.
      In correlation analysis, we found genes (LOC104968484, SLC25A6, GADD45G, BAX, XIAP, APAF1, ATM, MDM4, BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) that had a strong correlation with WBC count and W-SCC in SNT cows only. The relationship between these genes and blood routine test parameters suggested that the SNT phenotype affected the immunity of cows. The LOC104968484, SLC25A6, GADD45G, and BAX genes were negatively correlated with WBC count and W-SCC only in SNT cows, and the expression levels of these genes were significantly upregulated in SNT cows. Immunoglobulin heavy variable 4-59 (LOC104968484) and solute carrier family 25 member 6 (SLC25A6) were enriched in the ENT formation pathway. Growth arrest and DNA damage 45G (GADD45G) and B-cell lymphoma-2-associated X protein (BAX) promote apoptosis (
      • Edlich F.
      BCL-2 proteins and apoptosis: Recent insights and unknowns.
      ;
      • Shin G.-T.
      • Lee H.J.
      • Park J.E.
      Growth arrest and DNA damage 45 gamma is required for caspase-dependent renal tubular cell apoptosis.
      ). Apoptosis has been described as a potential mechanism involved in tissue destruction, which is caused by infectious and inflammatory diseases, such as mastitis (
      • Mi S.
      • Tang Y.
      • Shi L.
      • Liu X.
      • Si J.
      • Yao Y.
      • Augustino S.M.A.
      • Fang L.
      • Yu Y.
      Protective roles of folic acid in the responses of bovine mammary epithelial cells to different virulent Staphylococcus aureus strains.
      ).
      • Raj S.R.
      • Das D.N.
      • Mondal S.
      • Ashokan M.
      • Thota L.N.
      • Karuthadurai T.
      • Tej N.K.J.
      • Ramesha K.P.
      Expression analysis of pro-apoptotic BAX and anti-apoptotic BCL-2 genes in relation to lactation performance in Deoni and Holstein Friesian crossbred cows.
      showed that the BAX gene has an adverse effect on milk yield (
      • Raj S.R.
      • Das D.N.
      • Mondal S.
      • Ashokan M.
      • Thota L.N.
      • Karuthadurai T.
      • Tej N.K.J.
      • Ramesha K.P.
      Expression analysis of pro-apoptotic BAX and anti-apoptotic BCL-2 genes in relation to lactation performance in Deoni and Holstein Friesian crossbred cows.
      ). After heat stress, the relative survival rate of mammary epithelial cells decreases significantly, and the expression level of BAX increases significantly (
      • Yang M.
      • Wen H.
      • Mi S.
      • Yu Y.
      The influence of supernumerary teats on milk production performance and blood routine in Chinese Holsteins.
      ). Therefore, the upregulation of BAX in dairy cows may predict adverse results.
      Additionally, 4 genes (XIAP, APAF1, ATM, and MDM4) were positively correlated with WBC count and W-SCC only in SNT cows, and these genes were downregulated in SNT cows. X-linked inhibitor of apoptosis protein (XIAP) is a crucial anti-apoptotic protein, which is the only inhibitor of the apoptosis protein to directly interact with caspases (
      • Shang J.
      • Chen W.-M.
      • Wang Z.-H.
      • Wei T.-N.
      • Chen Z.-Z.
      • Wu W.-B.
      CircPAN3 mediates drug resistance in acute myeloid leukemia through the miR-153-5p/miR-183-5p-XIAP axis.
      ). Studies have shown that deficiency of the XIAP protein in humans leads to immunodeficiency and inflammatory bowel disease (
      • Rigaud S.
      • Fondanèche M.-C.
      • Lambert N.
      • Pasquier B.
      • Mateo V.
      • Soulas P.
      • Galicier L.
      • Le Deist F.
      • Rieux-Laucat F.
      • Revy P.
      • Fischer A.
      • de Saint Basile G.
      • Latour S.
      XIAP deficiency in humans causes an X-linked lymphoproliferative syndrome.
      ). Activation of T cells plays an important role in the clearance of foreign antigens. When antigens are removed, the activated T cells should be removed to avoid unnecessary inflammation. Apoptosomes in T cells play an important role in this process (
      • Reiner S.L.
      Decision making during the conception and career of CD4+ T cells.
      ), and apoptotic peptidase activating factor 1 (APAF1) participates in the formation and activation of the apoptosome (
      • Zou H.
      • Henzel W.J.
      • Liu X.
      • Lutschg A.
      • Wang X.
      Apaf-1, a human protein homologous to C. elegans CED-4, participates in cytochrome c-dependent activation of caspase-3.
      ). However, APAF1 was downregulated in SNT cows, which indicated that SNT cows were likely to suffer from unnecessary inflammation. ATM affects mammary duct development and mammary epithelial hyperplasia (
      • Bowen T.J.
      • Yakushiji H.
      • Montagna C.
      • Jain S.
      • Ried T.
      • Wynshaw-Boris A.
      Atm heterozygosity cooperates with loss of Brca1 to increase the severity of mammary gland cancer and reduce ductal branching.
      ). Mice without ATM failed to develop mature mammary glands (
      • Bowen T.J.
      • Yakushiji H.
      • Montagna C.
      • Jain S.
      • Ried T.
      • Wynshaw-Boris A.
      Atm heterozygosity cooperates with loss of Brca1 to increase the severity of mammary gland cancer and reduce ductal branching.
      ;
      • Dyer L.M.
      • Kepple J.D.
      • Ai L.
      • Kim W.-J.
      • Stanton V.L.
      • Reinhard M.K.
      • Backman L.R.F.
      • Streitfeld W.S.
      • Babu N.R.
      • Treiber N.
      • Scharffetter-Kochanek K.
      • McKinnon P.J.
      • Brown K.D.
      ATM is required for SOD2 expression and homeostasis within the mammary gland.
      ), indicating the importance of ATM to mammary gland development. MDM4 is an inhibitor of p53 (
      • Kawai H.
      • Lopez-Pajares V.
      • Kim M.M.
      • Wiederschain D.
      • Yuan Z.-M.
      RING domain-mediated interaction is a requirement for MDM2's E3 ligase activity.
      ), and the p53 protein induces apoptosis (
      • Wawryk-Gawda E.
      • Chylińska-Wrzos P.
      • Lis-Sochocka M.
      • Chłapek K.
      • Bulak K.
      • Jędrych M.
      • Jodłowska-Jędrych B.
      P53 protein in proliferation, repair and apoptosis of cells.
      ). These results suggested that SNT may affect the apoptosis and immunity of dairy cows. Numerous studies have shown that apoptosis, inflammation, and immunity exist almost simultaneously in most cases (
      • Chen Z.
      • Liang Y.-S.
      • Zong W.-C.
      • Guo J.
      • Zhou J.-H.
      • Mao Y.-J.
      • Ji D.-J.
      • Jiao P.-X.
      • Juan J.L.
      • Yang Z.-P.
      NH4Cl promotes apoptosis and inflammation in bovine mammary epithelial cells via the circ02771/miR-194b/TGIF1 axis.
      ;
      • Ma Y.
      • Ma X.
      • An Y.
      • Sun Y.
      • Dou W.
      • Li M.
      • Bao H.
      • Zhang C.
      Green tea polyphenols alleviate hydrogen peroxide-induced oxidative stress, inflammation, and apoptosis in bovine mammary epithelial cells by activating ERK1/2–NFE2L2-HMOX1 pathways.
      ;
      • Meng M.
      • Wang L.
      • Wang Y.
      • Ma N.
      • Xie W.
      • Chang G.
      • Shen X.
      A high-concentrate diet provokes inflammation, endoplasmic reticulum stress, and apoptosis in mammary tissue of dairy cows through the upregulation of STIM1/ORAI1.
      ).
      In PPI network analysis, we found 6 hub genes (BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) are also associated with SNT, and a close relationship was observed between the 6 hub genes and the IGR-SNT (LOC104968484, SLC25A6, GADD45G, BAX, XIAP, APAF1, ATM, and MDM4). The 6 hub genes were significantly, or tended to be significantly, positively correlated with WBC count and W-SCC only in SNT cows. BDP1, the hub gene with the highest score, is a subunit of transcription initiation factor IIIB (TFIIIB). Initiation of gene transcription by RNA polymerase III (Pol III) requires the activity of TFIIIB (
      • Gouge J.
      • Guthertz N.
      • Kramm K.
      • Dergai O.
      • Abascal-Palacios G.
      • Satia K.
      • Cousin P.
      • Hernandez N.
      • Grohmann D.
      • Vannini A.
      Molecular mechanisms of Bdp1 in TFIIIB assembly and RNA polymerase III transcription initiation.
      ). In addition to BDP1, CEP350 (
      • Kumar A.
      • Rajendran V.
      • Sethumadhavan R.
      • Purohit R.
      CEP proteins: The knights of centrosome dynasty.
      ) and MED13 (
      • Zhou W.
      • Cai H.
      • Li J.
      • Xu H.
      • Wang X.
      • Men H.
      • Zheng Y.
      • Cai L.
      Potential roles of mediator complex subunit 13 in cardiac diseases.
      ) contribute to the cell cycle. The TAOK1 protein negatively regulates IL-17-mediated signaling and inflammation (
      • Zhang Z.
      • Tang Z.
      • Ma X.
      • Sun K.
      • Fan L.
      • Fang J.
      • Pan J.
      • Wang X.
      • An H.
      • Zhou J.
      TAOK1 negatively regulates IL-17-mediated signaling and inflammation.
      ). The SMG1 protein mediates post-transcriptional regulation of genes (
      • Lloyd J.P.
      • Davies B.
      SMG1 is an ancient nonsense-mediated mRNA decay effector.
      ;
      • Roberts T.L.
      • Ho U.
      • Luff J.
      • Lee C.S.
      • Apte S.H.
      • MacDonald K.P.A.
      • Raggat L.J.
      • Pettit A.R.
      • Morrow C.A.
      • Waters M.J.
      • Chen P.
      • Woods R.G.
      • Thomas G.P.
      • St. Pierre L.
      • Farah C.S.
      • Clarke R.A.
      • Brown J.A.L.
      • Lavin M.F.
      Smg1 haploinsufficiency predisposes to tumor formation and inflammation.
      ). PIF1 is involved in DNA replication (
      • Blasiak J.
      • Szczepańska J.
      • Sobczuk A.
      • Fila M.
      • Pawlowska E.
      RIF1 links replication timing with fork reactivation and DNA double-strand break repair.
      ). The 6 hub genes that were downregulated in SNT cows were positively correlated with downregulated IGR-SNT and negatively correlated with upregulated IGR-SNT. Thus, a strong regulatory relationship exists between the hub genes of SG and IGR-SNT. The expression levels of IGR-SNT and the 6 hub genes all have a correlation with WBC count and W-SCC only in SNT cows. Therefore, these genes are associated with the SNT phenotype. However, we do not know the causal relationship between these genes and the SNT phenotype. That is, the presented link between these genes and SNT may be due to a higher infection load in SNT cows (SNT holding the bacteria). It is also possible that these genes are the cause of the SNT phenotype, but we did not conduct trials to prove this in the current study, which is a limitation of this study.
      Alternative splicing may mediate the effect of SNT on dairy cows and the formation of SNT. Alternative splicing, which plays an active role in normal growth and development of higher eukaryotes, contributes to cell differentiation, organ development, and tissue homeostasis (
      • Wang E.T.
      • Sandberg R.
      • Luo S.
      • Khrebtukova I.
      • Zhang L.
      • Mayr C.
      • Kingsmore S.F.
      • Schroth G.P.
      • Burge C.B.
      Alternative isoform regulation in human tissue transcriptomes.
      ). In the alternative splicing analysis, about one-fifth of SG were DSG, and 2 IGR-SNT (BAX and MDM4) were DSG. Thus, alternative splicing may mediate the effect of SNT on dairy cows or formation of SNT. The KEGG enrichment analysis of DSG also showed multiple inflammatory and immune-related pathways that were significantly enriched, such as the T cell receptor pathway, NOD-like receptor signaling pathway, MAPK signaling pathway, and chemokine signaling pathway. The mechanism of alternative splicing participating in the effect of SNT on dairy cows needs to be further studied.

      CONCLUSIONS

      Our study is the first to demonstrate that the SNT phenotype is positively associated with inflammation-related pathways, and negatively associated with the prolactin signaling pathway and carbohydrate digestion and absorption from the transcriptional level, by comparing blood transcriptome sequencing data of SNT and non-SNT cows. We screened genes (LOC104968484, SLC25A6, GADD45G, BAX, APAF1, ATM, XIAP, MDM4, BDP1, CEP350, MED13, TAOK1, SMG1, and RIF1) associated with SNT through correlation analysis. The pattern of correlation between genes (IGR-SNT and the 6 hub genes) and WBC count and W-SCC showed that the SNT phenotype will reduce the immunity of cows. Therefore, the SNT phenotype has an adverse effect on dairy cows, including decreased immunity and increased inflammation, and may affect lactation. The genes we found can be used as a reference for further searching for candidate genes of the SNT phenotype.

      ACKNOWLEDGMENTS

      This work was funded by the National Key R&D Program of China (grant nos. 2021YFD1200903 and 2021YFD1200900, Beijing, China), the NSFC-PSF joint project (grant no. 31961143009, Beijing, China), the Beijing Dairy Industry Innovation Team (grant no. BAIC06, Beijing, China), the Beijing Natural Science Foundation (grant no. 6182021, Beijing, China), and the Earmarked Fund for Modern Agro-industry Technology Research System (grant no. CARS-36). The authors thank all the members of the Animal Molecular and Quantitative Genetics Laboratory in China Agricultural University (Beijing, China) for their assistance in the study design and data analysis. The authors have not stated any conflicts of interest.

      REFERENCES

        • Batajoo K.K.
        • Shaver R.D.
        Impact of nonfiber carbohydrate on intake, digestion, and milk-production by dairy-cows.
        J. Dairy Sci. 1994; 77 (8083418): 1580-1588
        • Blasiak J.
        • Szczepańska J.
        • Sobczuk A.
        • Fila M.
        • Pawlowska E.
        RIF1 links replication timing with fork reactivation and DNA double-strand break repair.
        Int. J. Mol. Sci. 2021; 22 (34768871)11440
        • Bowen T.J.
        • Yakushiji H.
        • Montagna C.
        • Jain S.
        • Ried T.
        • Wynshaw-Boris A.
        Atm heterozygosity cooperates with loss of Brca1 to increase the severity of mammary gland cancer and reduce ductal branching.
        Cancer Res. 2005; 65 (16204043): 8736-8746
        • Brisken C.
        • Rajaram R.D.
        Alveolar and lactogenic differentiation.
        J. Mammary Gland Biol. Neoplasia. 2006; 11 (17111223): 239-248
        • Brka M.
        • Reinsch N.
        • Junge W.
        • Kalm E.
        Frequency and heritability of supernumerary teats in German Holsteins (Germany).
        Zuchtungskunde. 2000; 72: 17-27
        • Brka M.
        • Reinsch N.
        • Kalm E.
        Frequency and heritability of supernumerary teats in German Simmental and German Brown Swiss cows.
        J. Dairy Sci. 2002; 85 (12201539): 1881-1886
        • Brown J.
        • Schwartz R.A.
        Supernumerary nipples and renal malformations: A family study.
        J. Cutan. Med. Surg. 2004; 8 (15578130): 170-172
        • Butty A.M.
        • Frischknecht M.
        • Gredler B.
        • Neuenschwander S.
        • Moll J.
        • Bieber A.
        • Baes C.F.
        • Seefried F.R.
        Genetic and genomic analysis of hyperthelia in Brown Swiss cattle.
        J. Dairy Sci. 2017; 100 (27865493): 402-411
        • Cánovas A.
        • Rincón G.
        • Bevilacqua C.
        • Islas-Trejo A.
        • Brenaut P.
        • Hovey R.C.
        • Boutinaud M.
        • Morgenthaler C.
        • VanKlompenberg M.K.
        • Martin P.
        • Medrano J.F.
        Comparison of five different RNA sources to examine the lactating bovine mammary gland transcriptome using RNA-sequencing.
        Sci. Rep. 2014; 4 (25001089)5297
        • Chen Q.
        Supplemental tables and figures. figshare. Figure.
        • Chen Z.
        • Liang Y.-S.
        • Zong W.-C.
        • Guo J.
        • Zhou J.-H.
        • Mao Y.-J.
        • Ji D.-J.
        • Jiao P.-X.
        • Juan J.L.
        • Yang Z.-P.
        NH4Cl promotes apoptosis and inflammation in bovine mammary epithelial cells via the circ02771/miR-194b/TGIF1 axis.
        J. Integr. Agric. 2022; 21: 1161-1176
        • Davila J.I.
        • Fadra N.M.
        • Wang X.
        • McDonald A.M.
        • Nair A.A.
        • Crusan B.R.
        • Wu X.
        • Blommel J.H.
        • Jen J.
        • Rumilla K.M.
        • Jenkins R.B.
        • Aypar U.
        • Klee E.W.
        • Kipp B.R.
        • Halling K.C.
        Impact of RNA degradation on fusion detection by RNA-seq.
        BMC Genomics. 2016; 17 (27765019): 814
        • Dyer L.M.
        • Kepple J.D.
        • Ai L.
        • Kim W.-J.
        • Stanton V.L.
        • Reinhard M.K.
        • Backman L.R.F.
        • Streitfeld W.S.
        • Babu N.R.
        • Treiber N.
        • Scharffetter-Kochanek K.
        • McKinnon P.J.
        • Brown K.D.
        ATM is required for SOD2 expression and homeostasis within the mammary gland.
        Breast Cancer Res. Treat. 2017; 166 (28849346): 725-741
        • Edlich F.
        BCL-2 proteins and apoptosis: Recent insights and unknowns.
        Biochem. Biophys. Res. Commun. 2018; 500 (28676391): 26-34
        • Fraser D.
        Toward a global perspective on farm animal welfare.
        Appl. Anim. Behav. Sci. 2008; 113: 330-339
        • Gao X.
        • Oba M.
        Effect of increasing dietary nonfiber carbohydrate with starch, sucrose, or lactose on rumen fermentation and productivity of lactating dairy cows.
        J. Dairy Sci. 2016; 99 (26585468): 291-300
        • Gouge J.
        • Guthertz N.
        • Kramm K.
        • Dergai O.
        • Abascal-Palacios G.
        • Satia K.
        • Cousin P.
        • Hernandez N.
        • Grohmann D.
        • Vannini A.
        Molecular mechanisms of Bdp1 in TFIIIB assembly and RNA polymerase III transcription initiation.
        Nat. Commun. 2017; 8 (28743884): 130
        • Grayson P.C.
        • Kaplan M.J.
        At the bench: Neutrophil extracellular traps (NETs) highlight novel aspects of innate immune system involvement in autoimmune diseases.
        J. Leukoc. Biol. 2016; 99 (26432901): 253-264
        • Hersh J.
        Association of supernumerary nipples and renal anomalies—Reply.
        Arch. Pediatr. Adolesc. Med. 1988; 142 (3285662): 591-592
        • Hurley W.L.
        Review: Mammary gland development in swine: Embryo to early lactation.
        Animal. 2019; 13 (31280748): s11-s19
        • Joerg H.
        • Meili C.
        • Ruprecht O.
        • Bangerter E.
        • Burren A.
        • Bigler A.
        A genome-wide association study reveals a QTL influencing caudal supernumerary teats in Holstein cattle.
        Anim. Genet. 2014; 45 (25204440): 871-873
        • Kawai H.
        • Lopez-Pajares V.
        • Kim M.M.
        • Wiederschain D.
        • Yuan Z.-M.
        RING domain-mediated interaction is a requirement for MDM2's E3 ligase activity.
        Cancer Res. 2007; 67 (17616658): 6026-6030
        • Klopf J.
        • Brostjan C.
        • Eilenberg W.
        • Neumayer C.
        Neutrophil extracellular traps and their implications in cardiovascular and inflammatory disease.
        Int. J. Mol. Sci. 2021; 22 (33429925): 559
        • Kumar A.
        • Rajendran V.
        • Sethumadhavan R.
        • Purohit R.
        CEP proteins: The knights of centrosome dynasty.
        Protoplasma. 2013; 250 (23456457): 965-983
        • Lacasse P.
        • Lollivier V.
        • Bruckmaier R.M.
        • Boisclair Y.R.
        • Wagner G.F.
        • Boutinaud M.
        Effect of the prolactin-release inhibitor quinagolide on lactating dairy cows.
        J. Dairy Sci. 2011; 94 (21338795): 1302-1309
        • Liu B.
        • Li Q.
        • Gong Z.
        • Zhao J.
        • Gu B.
        • Feng S.
        Staphylococcus aureus lipoproteins play crucial roles in inducing inflammatory responses and bacterial internalization into bovine mammary epithelial cells.
        Microb. Pathog. 2022; 162 (34921958)105364
        • Lloyd J.P.
        • Davies B.
        SMG1 is an ancient nonsense-mediated mRNA decay effector.
        Plant J. 2013; 76 (24103012): 800-810
        • Lollivier V.
        • Lacasse P.
        • Angulo Arizala J.
        • Lamberton P.
        • Wiart S.
        • Portanguen J.
        • Bruckmaier R.
        • Boutinaud M.
        In vivo inhibition followed by exogenous supplementation demonstrates galactopoietic effects of prolactin on mammary tissue and milk production in dairy cows.
        J. Dairy Sci. 2015; 98 (26387019): 8775-8787
        • Ma Y.
        • Ma X.
        • An Y.
        • Sun Y.
        • Dou W.
        • Li M.
        • Bao H.
        • Zhang C.
        Green tea polyphenols alleviate hydrogen peroxide-induced oxidative stress, inflammation, and apoptosis in bovine mammary epithelial cells by activating ERK1/2–NFE2L2-HMOX1 pathways.
        Front. Vet. Sci. 2022; 8 (35146014)804241
        • Mackenzie K.J.J.
        • Marshall F.H.A.
        On the presence of supernumerary mammary glands in cows and on their functional activity.
        J. Agric. Sci. 1925; 15: 30-35
        • Martin P.
        • Palhiere I.
        • Tosser-Klopp G.
        • Rupp R.
        Heritability and genome-wide association mapping for supernumerary teats in French Alpine and Saanen dairy goats.
        J. Dairy Sci. 2016; 99 (27544860): 8891-8900
        • Méhes K.
        Association of supernumerary nipples with other anomalies.
        J. Pediatr. 1979; 95 (448569): 268-272
        • Meng M.
        • Wang L.
        • Wang Y.
        • Ma N.
        • Xie W.
        • Chang G.
        • Shen X.
        A high-concentrate diet provokes inflammation, endoplasmic reticulum stress, and apoptosis in mammary tissue of dairy cows through the upregulation of STIM1/ORAI1.
        J. Dairy Sci. 2022; 105 (35094865): 3416-3429
        • Mi S.
        • Tang Y.
        • Shi L.
        • Liu X.
        • Si J.
        • Yao Y.
        • Augustino S.M.A.
        • Fang L.
        • Yu Y.
        Protective roles of folic acid in the responses of bovine mammary epithelial cells to different virulent Staphylococcus aureus strains.
        Biology (Basel). 2021; 10 (34827157)1164
        • Palacios C.
        • Abecia J.-A.
        Supernumerary teat removal can be avoided in dairy sheep.
        J. Appl. Anim. Welf. Sci. 2014; 17 (24665954): 178-182
        • Pausch H.
        • Jung S.
        • Edel C.
        • Emmerling R.
        • Krogmeier D.
        • Götz K.-U.
        • Fries R.
        Genome-wide association study uncovers four QTL predisposing to supernumerary teats in cattle.
        Anim. Genet. 2012; 43 (22497297): 689-695
        • Peng W.F.
        • Xu S.S.
        • Ren X.
        • Lv F.H.
        • Xie X.L.
        • Zhao Y.X.
        • Zhang M.
        • Shen Z.Q.
        • Ren Y.L.
        • Gao L.
        • Shen M.
        • Kantanen J.
        • Li M.H.
        A genome-wide association study reveals candidate genes for the supernumerary nipple phenotype in sheep (Ovis aries).
        Anim. Genet. 2017; 48 (28703336): 570-579
        • Raj S.R.
        • Das D.N.
        • Mondal S.
        • Ashokan M.
        • Thota L.N.
        • Karuthadurai T.
        • Tej N.K.J.
        • Ramesha K.P.
        Expression analysis of pro-apoptotic BAX and anti-apoptotic BCL-2 genes in relation to lactation performance in Deoni and Holstein Friesian crossbred cows.
        Anim. Biotechnol. 2022; 22 (35067189): 1-8
        • Reiner S.L.
        Decision making during the conception and career of CD4+ T cells.
        Nat. Rev. Immunol. 2009; 9 (19172726): 81-82
        • Rigaud S.
        • Fondanèche M.-C.
        • Lambert N.
        • Pasquier B.
        • Mateo V.
        • Soulas P.
        • Galicier L.
        • Le Deist F.
        • Rieux-Laucat F.
        • Revy P.
        • Fischer A.
        • de Saint Basile G.
        • Latour S.
        XIAP deficiency in humans causes an X-linked lymphoproliferative syndrome.
        Nature. 2006; 444 (17080092): 110-114
        • Roberts T.L.
        • Ho U.
        • Luff J.
        • Lee C.S.
        • Apte S.H.
        • MacDonald K.P.A.
        • Raggat L.J.
        • Pettit A.R.
        • Morrow C.A.
        • Waters M.J.
        • Chen P.
        • Woods R.G.
        • Thomas G.P.
        • St. Pierre L.
        • Farah C.S.
        • Clarke R.A.
        • Brown J.A.L.
        • Lavin M.F.
        Smg1 haploinsufficiency predisposes to tumor formation and inflammation.
        Proc. Natl. Acad. Sci. USA. 2013; 110 (23277562): E285-E294
        • Gallego Romero I.
        • Pai A.A.
        • Tung J.
        • Gilad Y.
        RNA-seq: Impact of RNA degradation on transcript quantification.
        BMC Biol. 2014; 12 (24885439): 42
        • Shang J.
        • Chen W.-M.
        • Wang Z.-H.
        • Wei T.-N.
        • Chen Z.-Z.
        • Wu W.-B.
        CircPAN3 mediates drug resistance in acute myeloid leukemia through the miR-153-5p/miR-183-5p-XIAP axis.
        Exp. Hematol. 2019; 70 (30395908): 42-54.e3
        • Shin G.-T.
        • Lee H.J.
        • Park J.E.
        Growth arrest and DNA damage 45 gamma is required for caspase-dependent renal tubular cell apoptosis.
        PLoS One. 2019; 14 (30794682)e0212818
        • Skjervold H.
        Supernumerary teats in cattle.
        Hereditas. 1960; 46: 71-74
        • Vasileiou N.G.C.
        • Mavrogianni V.S.
        • Petinaki E.
        • Fthenakis G.C.
        Predisposing factors for bacterial mastitis in ewes.
        Reprod. Domest. Anim. 2019; 54 (31361921): 1424-1431
        • Villanueva E.
        • Yalavarthi S.
        • Berthier C.C.
        • Hodgin J.B.
        • Khandpur R.
        • Lin A.M.
        • Rubin C.J.
        • Zhao W.
        • Olsen S.H.
        • Klinker M.
        • Shealy D.
        • Denny M.F.
        • Plumas J.
        • Chaperot L.
        • Kretzler M.
        • Bruce A.T.
        • Kaplan M.J.
        Netting neutrophils induce endothelial damage, infiltrate tissues, and expose immunostimulatory molecules in systemic lupus erythematosus.
        J. Immunol. 2011; 187 (21613614): 538-552
        • Wang E.T.
        • Sandberg R.
        • Luo S.
        • Khrebtukova I.
        • Zhang L.
        • Mayr C.
        • Kingsmore S.F.
        • Schroth G.P.
        • Burge C.B.
        Alternative isoform regulation in human tissue transcriptomes.
        Nature. 2008; 456 (18978772): 470-476
        • Wawryk-Gawda E.
        • Chylińska-Wrzos P.
        • Lis-Sochocka M.
        • Chłapek K.
        • Bulak K.
        • Jędrych M.
        • Jodłowska-Jędrych B.
        P53 protein in proliferation, repair and apoptosis of cells.
        Protoplasma. 2014; 251 (24043441): 525-533
        • Wen H.
        • Luo H.
        • Mi S.
        • Liu X.
        • Wang Y.
        • Xiao W.
        • Yu Y.
        Effect of supernumerary teat on milk performance and genetic expression analysis of LGR5 gene in dairy cows.
        China Animal Husbandry & Veterinary Medicine. 2019; 46: 3642-3649
        • Wen H.
        • Luo H.
        • Yang M.
        • Augustino S.M.A.
        • Wang D.
        • Mi S.
        • Guo Y.
        • Zhang Y.
        • Xiao W.
        • Wang Y.
        • Yu Y.
        Genetic parameters and weighted single-step genome-wide association study for supernumerary teats in Holstein cattle.
        J. Dairy Sci. 2021; 104 (34482976): 11867-11877
        • Yang M.
        • Wen H.
        • Mi S.
        • Yu Y.
        The influence of supernumerary teats on milk production performance and blood routine in Chinese Holsteins.
        Zhongguo Xumu Zazhi. 2022; 58: 69-74
        • Zhang L.
        • Peng F.
        • Yu F.
        • Wan L.
        • Zhou Z.Q.
        Expression of ESR1, PRLR, GHR, and IGF1R in mammary glands of Hu sheep with four teats.
        Czech J. Anim. Sci. 2019; 64: 49-58
        • Zhang Z.
        • Tang Z.
        • Ma X.
        • Sun K.
        • Fan L.
        • Fang J.
        • Pan J.
        • Wang X.
        • An H.
        • Zhou J.
        TAOK1 negatively regulates IL-17-mediated signaling and inflammation.
        Cell. Mol. Immunol. 2018; 15 (29400705): 794-802
        • Zhang Z.
        • Wei Q.
        • Zeng Y.
        • Jia X.
        • Su H.
        • Lin W.
        • Xing N.
        • Bai H.
        • He Y.
        • Wang Q.
        Effect of Hordei Fructus Germinatus on differential gene expression in the prolactin signaling pathway in the mammary gland of lactating rats.
        J. Ethnopharmacol. 2021; 268 (33217517)113589
        • Zhou W.
        • Cai H.
        • Li J.
        • Xu H.
        • Wang X.
        • Men H.
        • Zheng Y.
        • Cai L.
        Potential roles of mediator complex subunit 13 in cardiac diseases.
        Int. J. Biol. Sci. 2021; 17 (33390853): 328-338
        • Zou H.
        • Henzel W.J.
        • Liu X.
        • Lutschg A.
        • Wang X.
        Apaf-1, a human protein homologous to C. elegans CED-4, participates in cytochrome c-dependent activation of caspase-3.
        Cell. 1997; 90 (9267021): 405-413