Functional Analysis of Haplotypes in Bovine PSAP Gene and Their Relationship with Beef Cattle Production Traits

Simple Summary With the rapid development of information technology and molecular biotechnology, animal molecular breeding technology is playing an increasingly important role in beef cattle breeding. Prosaposin (PSAP) is involved in regulating the growth and development of animals, and it is reported that PSAP is an important marker-assisted selection (MAS) in cattle herd. The purpose of this study was to explore the novel variants in 3’ UTR of cattle PSAP and evaluate their effects on the morphological traits of four Chinese cattle breeds. In this study, 13 variants were identified in the PSAP 3’ UTR from 501 individuals belonging to four cattle breeds. In Nanyang cattle, the distribution of haplotypes was different from the other three breeds. Two groups of haplotypes had association with morphological traits by changing the secondary structures of PSAP 3’ UTR rather than the miR-184 target sites. This study not only expands the genetic variation spectrum of cattle PSAP but also contributes to MAS genetics and breeding of Chinese cattle breeds. Abstract The purpose of this study was to explore functional variants in the prosaposin (PSAP) three prime untranslated region (3’ UTR) and clarify the relationship between the variants and morphological traits. Through Sanger sequencing, 13 variations were identified in bovine PSAP in four Chinese cattle breeds, with six of them being loci in 3’ UTR. In particular, Nanyang (NY) cattle had a special genotype and haplotype distribution compared to the other three breeds. NY cattle with ACATG and GCGTG haplotypes had higher morphological traits than GTACA and GTACG haplotypes. The results of dual-luciferase reporter assay showed that ACATG and GCGTG haplotypes affected the morphological traits of NY cattle by altering the secondary structure of PSAP 3’ UTR rather than the miR-184 target sites. The findings of this study could be an evidence of a complex and varying mechanism between variants and animal morphological traits and could be used to complement candidate genes for molecular breeding.


Introduction
With continuous economic development and growth of formalized cultivation, beef cattle breeding is in urgent need of strengthening. Traditional breeding methods are slow to produce results, making it difficult to meet demand. Fortunately, over the last decade, bioinformatics has revolutionized livestock breeding. Compared with traditional breeding, molecular breeding has a number of advantages, such as saving time and shortening the generation interval [1,2]. For example, genome-wide association studies (GWAS) have mapped thousands of genetic variants associated with animal development, which is an unprecedented high-resolution genetic characterization of animal breeding [3][4][5][6]. However, growing evidence suggests that successful livestock breeding requires deep understanding of the regulatory mechanism between the genotype and phenotype relationship [7][8][9][10]. The location of these variants can result in phenotype differences, including the level of DNA

Materials and Methods
All experimental procedures were performed in accordance with the Regulations for the Administration of Affairs Concerning Experimental Animals approved by the State Council of People's Republic of China. The study was approved by The Institutional Animal Care and Use Committee of Northwest A&F University (Permit Number: NWAFAC1019).

Samples and Data Collection
In order to explore variations in bovine PSAP 3' UTR, 501 cattle samples from four breeds (LX, n = 104, female; QC, n = 123, female; NY, n = 137, female; JX, n = 137, female) were collected (Table 1). All animals within a breed were managed in the same condition, and sufficient feed was provided by total metabolic rate (TMR). Their morphological traits were measured, including body weight (BW), body height (BH), body length (BL), chest circumference (ChC), chest depth (ChD), chest width (ChW), hucklebone width (HuW), and hip width (HW). Genome DNA was extracted from the leukocytes of venous blood by the phenol chloroform method [16,20]. All DNA samples underwent quality assurance and were uniformly diluted to 50 ng/µL and stored at −80 • C [16].

Primer Design and Variation Genotyping
According to the cattle reference genome (ARS-UCD1.2) from the NCBI database (http://www.ncbi.nlm.nih.gov/), a pair of specific primers was designed by NCBI primer blast for PCR amplification of cattle PSAP 3' UTR ( Table 2). The forward primer was located in the 14th intron, and the reverse primer was located in the downstream region of PSAP. PCR amplification was performed in 20 µL final volume containing 10 µL 2 × PCR mix, 0.5 µM forward and reverse primers, 50 ng genomic DNA, and ddH 2 O up to 20 µL. The PCR protocol was touchdown PCR as follows: 5 min predegeneration, followed by 10 cycles at 95 • C for 30 s, 60 • C for 30 s (starting at 60 • C and decreasing by 1 • C per cycle), 72 • C for 30 s, 25 cycles at 95 • C for 30 s, annealing temperature for 30 s, 72 • C for 30 s, and finally extending for 10 min. The PCR products of all individuals were sequenced by Sanger sequencing (Sangon, Shanghai, China) to determine the variation types.

Dual-Luciferase Reporter Assay
Human embryonic kidney cell line (HEK293T) was cultured in Dulbecco's modified Eagle medium (DMEM; Hyclone, Logan, UT, USA) supplemented with 10% fetal bovine serum (BI, Kibbutz Beit Haemek, Israel), 500 U/mL penicillin, and 100 µg/mL streptomycin (Hyclone, Logan, UT, USA). The cell incubator was controlled at 37 • C and 5% of the CO 2 volume. Four PSAP 3' UTR haplotypes were cloned into psiCHECK2 plasmid by complementary double-stranded annealing, vector double enzyme restriction, and adapter ligation ( Table 1). The cells were cotransfected with a mixture of 500 ng recombined reporter vectors and 3 µL miR-184 mimics. After 48 h, the luciferase activity was measured with a dual-luciferase reporter assay system (Promega, Madison, WI, USA). In the luciferase assay, the mimic NC was the negative control. The fold change of luciferase was calculated by comparing each miRNA to NC according to the manufacturer's instructions. Primers (Table 2) were synthesized by Sangon Biotech (Shanghai, China) Co., Ltd.

Statistical Analysis
Sequences were contrasted and analyzed by Snapgene (GSL Biotech, Chicago, IL, USA). All the population genetic data were calculated and analyzed using the website www.Msrcall.com, including the Hardy-Weinberg equilibrium (HWE), homozygosity (Ho), heterozygosity (He), effective allele number (Ne), and polymorphism information content (PIC). Every cattle breed was collected in the same farm and were of the same gender (female), had a similar age (2-6 years old), had the same feeding management, and had no genetic relationship within a population. The records of morphological traits were analyzed in each cattle breed independently. Mixed linear model analysis was used to establish the influence of different parameters on morphological traits, not including the effects of farm, sex, or age of dam and sire, which had no significant effects on the variation of traits in the four cattle populations in this study. The least squares mean was utilized for morphological traits among the different genotypes and haplotypes: Yi is the phenotypic value of morphological traits; µ is the overall population mean; Gj is the genotype; and ei is the random error. The chi-square test and ANOVA were used with SPSS software 18.0 (IBM, San Francisco, CA, USA) to test the genotype distribution and the association with morphological traits in different cattle breeds. The least significant difference (LSD) was used for multiple comparison. Statistical significance was measured at p < 0.05 and p < 0.01. Haplotype and linkage disequilibrium (LD) analyses were done using the website http://analysis.bio-x.cn/myAnalysis.php, and LD heatmap (R) and mRNA secondary structure were predicted using the online software RNAfold.

Population Parameters of 13 Variations in Four Cattle Breeds
The population parameters of four cattle breeds were calculated, including genotype frequency, allele frequency, HWE, Ho, He, Ne, and PIC (Tables 3 and 4). For the P3 locus, CC was the dominant genotype in LX, NY, and JX cattle, but TT was the dominant genotype in QC cattle (Table 3). Chi-square test ( Figure 2) showed the difference in allele and genotype frequency between QC cattle and the other three breeds more clearly. In the P8 locus, TT was the dominant genotype in NY, while CC was the dominant genotype in the other three groups (Table 4). Figure 2 more clearly shows the difference in genotype and allele frequency between NY and the other three breeds at locus P8. Incidentally, for the N6 (P7) locus, the quantity of TT and TC genotypes were not enough for statistical calculation in the present study ( Table 4). The rest of the genotype distribution of the 11 loci were the same in the four cattle breeds (Tables 3 and 4). Linkage disequilibrium analysis showed that the linkage of 13 loci in QC, LX, and JX cattle was different from the NY cattle ( Figure 3).

Association between 13 Variations and Morphological Traits
The associations between the 13 novel variants and cattle morphological traits were investigated. Seven of them were associated with morphological traits: P1, P3, P4, N1, N12, N31, and P12. For the P1 locus, JX cattle with AA genotype had higher waist width and hucklebone width than GG genotype, but NY cattle with GG genotype had higher body length than AA genotype (p < 0.05). For the P3 locus, LX cattle with CC genotype had higher abdominal circumference than TT genotype, and the TC genotype was the best genotype (p < 0.05). For the P4 locus, NY cattle with AA genotype had higher waist width than GG genotype, but LX cattle with GG genotype had higher abdominal circumference than AA genotype, and the GA genotype was the best genotype (p < 0.05). For the N1 locus, NY cattle with GG and AA genotypes had higher abdominal circumference than GA genotype (p < 0.05). For the N12 locus, JX cattle with AA genotype had higher hucklebone width than GA genotype (p < 0.05). For the N31 locus, NY cattle with AA and GG genotypes had higher abdominal circumference than GA genotype, while JX cattle with AA and GG genotypes had higher body length, waist width, and hucklebone width than GA genotype (p < 0.05). For the P12 locus, NY cattle with GG genotype had higher abdominal circumference than CC genotype, and LX cattle with GG genotype had higher waist width and body weight than AA genotype (p < 0.05) ( Table 5).

Population Parameters of 13 Variations in Four Cattle Breeds
The population parameters of four cattle breeds were calculated, including genotype frequency, allele frequency, HWE, Ho, He, Ne, and PIC (Tables 3 and 4). For the P3 locus, CC was the dominant genotype in LX, NY, and JX cattle, but TT was the dominant genotype in QC cattle (Table 3). Chi-square test (Figure 2) showed the difference in allele and genotype frequency between QC cattle and the other three breeds more clearly. In the P8 locus, TT was the dominant genotype in NY, while CC was the dominant genotype in the other three groups (Table 4). Figure 2 more clearly shows the difference in genotype and allele frequency between NY and the other three breeds at locus P8. Incidentally, for the N6 (P7) locus, the quantity of TT and TC genotypes were not enough for statistical calculation in the present study ( Table 4). The rest of the genotype distribution of the 11 loci were the same in the four cattle breeds (Tables 3 and 4). Linkage disequilibrium analysis showed that the linkage of 13 loci in QC, LX, and JX cattle was different from the NY cattle ( Figure 3).       177.77 a ± 6.95 (n = 11) 0.009

Association of Four Haplotypes of PSAP 3' UTR and Morphological Traits in NY Cattle
Considering the low frequency of the N6 locus, haplotypes of PSAP 3' UTR in four cattle breeds were structured with five SNP loci (N1, N8, N12, N13, and N31). For the N8 locus, the number of allele T and genotype TT were different between NY cattle and the other cattle breeds. Based on the differences in genotype and allele, the haplotype distribution was different between NY cattle and the other cattle breeds. The dominant haplotypes of NY cattle were GTACG, GTACA, CCGTG, and ACATG, but the dominant haplotypes of JX, LX, and QC cattle were GCACG, GCACA, ATATG, and GTGTG ( Figure 4).
The relationship between the haplotypes and morphological traits were analyzed in the four cattle breeds. Interestingly, the association between haplotypes and morphological traits only appeared in the NY cattle population. NY cattle with haplotypes ACATG and GCGTG had higher body height, body length, chest circumference, and body weight than those with haplotypes GTACA and GTACG (Table 6).
3' UTR is an important regulatory region of mRNA stability, and miRNAs widely participate in the regulation process. Based on the miRBase database, bta-miR-184 was found to have different binding capacity with the four haplotypes in NY cattle. Four PSAP haplotypes were cloned into psiCHECK2 plasmid. miR-184 mimics, mimics NC, and four haplotypes were cotransfected in the HEK293-T cell line. In the luciferase assay, the relative luciferase of the mimics NC group was higher than the miR-184 mimics group, and the miR-184 did not affect the targeting of PSAP 3' UTR. More interestingly, the relative luciferase between the two groups was different, with the relative luciferase of the ACATG and GCGTG group being higher than the GTACA and GTACG group ( Figure 5). Based on the evidence that SNPs could affect the mRNA secondary structure and protein folding, we concluded that two groups of haplotypes in PSAP 3' UTR could contribute to mRNA stability. Through a mRNA secondary structure prediction of four haplotypes by RNAfold, we found that the free energy of the ACATG and GCGTG group was higher than that of the GTACA and GTACG group ( Figure 6). This could be the reason for the differences in relative luciferase between the two groups. locus, the number of allele T and genotype TT were different between NY cattle other cattle breeds. Based on the differences in genotype and allele, the haplotyp bution was different between NY cattle and the other cattle breeds. The dominan types of NY cattle were GTACG, GTACA, CCGTG, and ACATG, but the domin lotypes of JX, LX, and QC cattle were GCACG, GCACA, ATATG, and GTGTG (F The relationship between the haplotypes and morphological traits were ana the four cattle breeds. Interestingly, the association between haplotypes and mor cal traits only appeared in the NY cattle population. NY cattle with haplotypes and GCGTG had higher body height, body length, chest circumference, and bod than those with haplotypes GTACA and GTACG (Table 6).   137.00 ± 1.96 a (n = 4) 129.58 ± 1.40 b (n = 12) 128.02 ± 0.59 b (n = 66) 0.001

Discussion
It is more strategically important than ever to preserve as much livestock diversity as possible to ensure enough basic gene reserve for future breeding needs [21]. There is less use of beef cattle genetic resources in the beef industry, and most of them are protected by the government [22,23]. Irrespective of the size of the farm, the large majority of breeds are commercial breeds rather than indigenous breeds [24]. Based on the complexity of China's geography, there is an extremely large amount of genetic resources that could be discovered and utilized from numerous cattle breeds. Molecular breeding is an appropriate choice, and analyzing inherent mechanism is necessary for their popularization and application [25].
PSAP has an important function in the Caspase, MAPK, PI3K/Akt, and TGF-β pathways, which have potential function in animal breeding [17]. For example, Guo found six SNPs in the 10th intron and 11th exon of PSAP that influenced carcass and meat quality traits [18]. In this study, 13 variations belonging to the 14th intron, 15th exon, 3' UTR, and downstream region of PSAP were found in four cattle breeds. After association analysis, we consider that the P4 locus (missense mutation) and six SNPs in 3' UTR could be seen as functional variants in beef breeding; however, there was not enough evidence to suggest this. The differences in genotype and allele distribution between NY cattle and the other three cattle breeds attracted our attention. The association between the haplotypes and body height, body length, chest circumference, and body weight were investigated, and we found the morphological traits with the ACATG and GCGTG group was higher than that with the GTACA and GTACG group. This could be regarded as an important basis for a potentially valuable quantitative trait locus (QTL).
The function of variants depends on their location. In the upstream, it can affect DNA methylation, the location of the CpG island, and transcription element recognition sites [26]. In the exon, it can affect protein translation and the spatial structure of proteins. In the intron, it can participate in the function of unknown coding genes, noncoding RNA, and gene alternative splicing. In the 3' UTR, it can affect targeting with miRNA. miRNA is a kind of small RNA that can target the gene 3' UTR and construct RNA-induced silencing complex (RISC) [27]. Many studies have demonstrated it as an important mechanism of variants. SNP-miRNA-mRNA interaction networks have an important role in human health and animal production [28,29]. Besides, all the variants in mRNA can affect stability due to the mRNA secondary structure and free energy [30,31]. However, the functional mechanisms underlying the associations are largely unknown [32][33][34]. To more clearly and comprehensively reveal the complex connections between the four haplotypes and miR-184, luciferase assay was performed to investigate the influence of the functional SNP loci. It is a pity that targeting miR-184 was not a mechanism of the functional SNPs in NY cattle PSAP 3' UTR. Compared with the expression of miRNA in muscle and adipose, its expression was not sufficient to mediate such a large effect on the morphological traits in NY cattle [34,35]. However, we found an unexpected phenomenon that the luciferase activity of the ACATG and GCGTG group was higher than the GTACA and GTACG group. To our surprise, the free energy of the ACATG and GCGTG group was lower than the GTACA and GTACG group. We hypothesize that the free energy of the haplotypes in the two groups is the reason for the differences in morphological traits in NY cattle.
More and more novel functional variants are being confirmed by the development of bioinformatics, but most of them lack underlying molecular mechanisms [36,37]. It is an important but difficult work to identify their potential value in animal breeding. Briefly, we identified 13 variations of loci in bovine PSAP, including six SNP loci in 3' UTR. The distribution of haplotypes in NY cattle was different from the other three breeds. Two groups of haplotypes had association with morphological traits by changing the secondary structures of PSAP 3' UTR rather than the miR-184 target sites. Our findings could be evidence of a complex and varying mechanism between variants and animal morphological traits and might be a useful MAS for beef cattle breeds.

Conclusions
In this study, 13 variants were identified in PSAP 3 UTR from 501 individuals of four Chinese cattle breeds, including three loci in the 14th intron, two loci in the 14th exon, six loci in the 3' UTR, and two loci in the downstream region of PSAP. NY cattle had a particular haplotype distribution compared to the other three cattle breeds. The results of the association analysis showed that NY cattle with haplotypes ACATG and GCGTG had better production traits than those with haplotypes GTACA and GTACG by changing the secondary structures of PSAP 3' UTR rather than the miR-184 target sites.

Informed Consent Statement: Not applicable.
Data Availability Statement: Data sharing not applicable. No new data were created or analyzed in this study. Data sharing is not applicable to this article.