Determination of the effect of functional single-nucleotide polymorphisms associated with glycerolipid synthesis on intramuscular fat deposition in Korean cattle steer

Abstract Intramuscular fat deposition in the longissimus dorsi muscle (LM) of Korean cattle steer is regulated by several genes related to lipid metabolism. One of these genes encodes the enzyme bovine glycerol-3-phosphate acyltransferase, mitochondrial (GPAM), which is located on the mitochondrial outer membrane and catalyzes the initial and committed step of glycerolipid synthesis in lipid metabolism of cattle. Previous studies have shown that the 3′-untranslated region (UTR) of the GPAM is quite extended and contains a polyadenylation signal site, erythroid 15-lipoxygenase differentiation control elements (15-LOX-DICEs), and cytoplasmic polyadenylation elements (CPEs) that affect the regulation of triacylglycerol synthesis. Therefore, the aim of this study was to identify single-nucleotide polymorphisms (SNPs) related to the regulation of glycerolipid synthesis in the 3′-UTR of GPAM and to verify the function of SNPs affecting the deposition of intramuscular fat in Korean cattle steer. In the present study, 11 SNPs were discovered in the 3′-UTR of GPAM. Among these SNPs, g.54853A>G, g.55441A>G, and g.55930C>T were significantly associated with marbling score in a Korean cattle steer population and were strongly correlated with each other within the GPAM gene. Furthermore, based on the results predicted by the RNAhybrid program, four putative microRNAs (miRNAs) were identified, and the above SNPs were found to present in the seed region of these miRNAs. These miRNAs have a differential binding affinity for each allele of SNPs g.54853A>G, g.55441A>G, and g.55930C>T. The in vivo evidence of intramuscular fat deposition in the LM tissue showed that these SNPs affected the regulation of intramuscular fat deposition in Korean cattle steer. Thus, the g.54853A>G, g.55441A>G, and g.55930C>T could be considered as causal mutations regulating intramuscular fat deposition in Korean cattle steer.


Introduction
Intramuscular fat deposition in the longissimus dorsi muscle (LM) of Korean cattle is regulated by several genes related to lipid metabolic processes, such as adipogenesis, lipogenesis, glycerolipid synthesis, and lipolysis. Several previous studies have reported that by increasing the deposition of intramuscular fat in Korean cattle steer, the expression levels of mRNAs related to adipogenesis, glycerolipid synthesis, and lipogenesis were upregulated, whereas those of mR-NAs related to lipolysis were downregulated (Jeong et al., 2012;Kim et al., 2008). Among these metabolic processes, the mRNA abundance of the glycerol-3-phosphate acyltransferase 1 (GPAT1), which catalyzes the initial and committed H. Kim et al.: Determination of the effect of functional SNPs step in glycerolipid biosynthesis, showed the greatest correlation with intramuscular fat content (Jeong et al., 2012).
Bovine glycerol-3-phosphate acyltransferase, mitochondrial (GPAM), also known as GPAT1, is a member of the GPAT gene family and is an enzyme that synthesizes lysophosphatidic acid (LPA) by transferring acyl groups to glycerol-3-phosphate (Yu et al., 2017). Thus, this enzyme catalyzes the initial and committed step in glycerolipid biosynthesis and plays a key role in regulating the level of cellular triacylglycerol in cattle (Yu et al., 2017).
The GPAM gene is located on bovine chromosome 26q22 and is composed of 21 exons and 20 introns. It is 3689 bp in length and has a much extended 3 -untranslated region (3 -UTR), which is longer than that of the other bovine genes. Furthermore, the 3 -UTR of GPAM contains a polyadenylation signal site, erythroid 15-lipoxygenase differentiation control elements (15-LOX-DICEs), and cytoplasmic polyadenylation elements (CPEs) that affect the regulation of triacylglycerol synthesis. Roy et al. (2006) reported that an extended 3 -UTR is significant for GPAM gene regulation, and a cis-element might be implicated in mRNA stabilization and translational control. The transmembrane domain of the GPAM protein is bound to the outer membrane of mitochondria and its N-and C-termini are located in the cytoplasm (Roy et al., 2006).
Recently, Yu et al. (2017) reported that the knockdown of GPAM expression significantly reduced the synthesis of triglycerides in bovine embryonic fibroblast (BEF) cells, and the genetic variation of GPAM was significantly associated with the fatty acid composition of intramuscular fat in cattle. Furthermore, Roy et al. (2006) reported that the cisregulatory elements (15-LOX-DICEs and CPEs), which regulate the genes related to lipid degradation and polyadenylation signal, respectively, are located in the 3 -UTR of the GPAM gene.
Therefore, the aim of this study was to identify SNPs related to the regulation of glycerolipid synthesis in the 3 -UTR of GPAM and to verify the function of SNPs affecting intramuscular fat deposition in Korean cattle steer.

Animals, DNA extraction, SNP discovery
Animal welfare issue was followed according to approved guidelines of the Animal Care and Use Committee of Hankyong National University. LM tissue samples were collected from Korean cattle (n = 192) raised in Pyeongchang (Gangwon, Republic of Korea). The marbling grade of this beef was classified according to the carcass grading standard of the Korea Institute of Animal Products Quality Evaluation (KAPE, 2017). The KAPE provided grades for beef marbling standard (BMS) score. All steers were maintained under constant environmental conditions, with two types of commercial feeds in six feedlots. Genomic DNA was extracted from the LM tissue using a LaboPass™ tissue mini kit (Cosmo Genetech, Seoul, Republic of Korea). In order to discover SNPs, the 3 -UTR sequence of GPAM was obtained from the National Center for Biotechnology Information (NCBI) GenBank database (accession no. NC_037353.1). The primers were designed using NCBI Primer-BLAST based on the selected polymorphism sites, and the primer sequences are shown in Table S1 in the Supplement. The sequencing was performed according to a previous study (Lee et al., 2010), and SNPs were discovered using the "SNP Hunting" option of the Sequencer v5.2.4 program (Gene Codes Corp., Ann Arbor, MI, USA). In order to map the functional SNPs on DNA, mRNA (NM_001012282.1), and protein (NP_01011282.1), these sequences were aligned using the NCBI graphical sequence viewer. The correlation coefficients between SNP pairs were calculated by Haploview 4.1 (Broad Institute, Cambridge, MA, USA) with genotypes of the GPAM gene in Korean cattle.

SNP genotyping and statistical analysis
SNPs were genotyped commercially using the Fluidigm ® SNP™-type assay platform according to a previous study (Oh et al., 2018). In order to evaluate the association SNPs and carcass traits, these data were analyzed using a generalized linear model (GLM) in SPSS v22 (IBM, Chicago, IL, USA) with the following equation: where Y ij kl is the observed marbling score of Korean cattle; µ is the overall mean; SNP k is the fixed effect of SNP genotype or haplotype; Farm i is the fixed effect of the feed type in farm; Sire j is the random effect of the sire; β age is the covariation of age; and e ij kl is random error. The correlation coefficient between SNP pairs in the 3 -UTR of GPAM was analyzed using the Haploview program (Broad Institute, USA). The relationship between SNPs and beef quality grades (1 ++ and 2) was analyzed using Fisher's exact test in SPSS v22 (IBM, USA).

Identification of SNPs and their functional characterization
In the present study, we discovered 11 SNPs within the 3 -UTR of GPAM in a Korean cattle steer population using the direct sequencing method. The positions of SNPs within the 3 -UTR of GPAM and their pairwise correlation coefficients are shown in Fig. 1. As shown in Fig. 1a, although the detected 11 SNPs were present in the 3 -UTR of GPAM, they were not specifically present in the sequences of four 15-LOX-DICEs and two CPEs present in the 3 -UTR. For determining linkage disequilibrium (LD) among 11 polymorphic SNPs in the 3 -UTR of GPAM, we calculated the correlation coefficient between SNP pairs using the Haploview software. The correlation coefficient between SNP pairs is shown in Fig. 1b. As shown in Fig. 1b, two LD blocks were detected in the 3 -UTR of GPAM. The LD block structure included two SNPs, whereas the other structure included four SNPs. Three SNPs, g.54853A>G, g.55441A>G, and g.55930C>T, were strongly correlated with each other, except for the g.55517G>A SNP.
Previous studies have reported that the 3 -UTR of the GPAM gene is longer than the open reading frame (ORF) region and includes a cis-element that plays a key role in mRNA stabilization and translational control (Roy et al., 2006;Yu et al., 2017). Furthermore, Yu et al. (2017) showed that the SNPs in the 3 -UTR of GPAM were significantly associated with fatty acid composition of intramuscular fat and marbling score in a beef cattle population.

Association of functional SNPs with marbling score
In order to evaluate the function of 11 SNPs in the regulation of glycerolipid synthesis, we analyzed the association between these SNPs and marbling score. The effects of these SNPs and their combinations on the marbling score of Korean cattle are shown in Tables 1 and 2.
As shown in Table 1, four SNPs (g.54316A>C, g.54853A>G, g.55441A>G, and g.55930C>T) were sig-   Block 2 GG-GG-TT 60 6.019 ± 0.269 0.044 GG-GG-TC 1 6.000 ± 1.749 GA-GA-TC 77 6.433 ± 0.234 AA-GA-CC 1 5.000 ± 1.749 AA-AA-TC 1 6.000 ± 1.749 AA-AA-CC 49 5.628 ± 0.295 nificantly associated with the marbling score of the Korean cattle steer population. The marbling scores in group with the heterozygote genotypes of SNPs g.54316A>C, g.54853A>G, g.55441A>G, and g.55930C>T were significantly higher than those with the homozygote genotypes of these SNPs. Furthermore, as shown in Table 2, the combination types of SNPs g.54853A>G, g.55441A>G, and g.55930C>T were significantly associated with the marbling score of the Korean cattle steer population. Especially, the group with the combination types GA, GA, and TC had the highest average marbling score in the studied Korean cattle steer population.

miRNA prediction and its target alteration by SNP alleles
In the present study, we predicted miRNAs that bind to candidate functional SNPs using a bioinformatics tool, in order to determine whether these SNPs had an effect on the regulation of glycerolipid synthesis. We identified four miRNAs that bind to candidate functional SNPs present in the seed region of miRNA using the online software programs Tar-getScan and miRNA_Target (Table 3). To increase the credibility of our results, we used the RNAhybrid software to quantitatively determine the binding affinity between miR-NAs and SNPs located in the seed region. As shown in Table 3, three SNPs (g.54853A>G, g.55441A>G, and g.55930C>T) were predicted to bind to four miRNAs (bta-miR-2418, bta-miR-375, bta-miR-2479, and bta-miR-2468) according to the seed region type. The seed region or seed sequence, which is two to seven nucleotides long at the 5 end of the miRNA sequence, is essential for the binding of the partially complementary miRNA to the mRNA. The seed region is classified into atypical sites, canonical sites and marginal depending on the number of nucleotides matching between the seed sequence and the mRNA (Bartel, 2009). Another important factor on which miRNA binding depends is whether is adenine or guanine is the first nucleotide at the 5 end of miRNA (Bartel, 2009). In addition, the binding efficiency is in the order of (a) 8mer, (b) 7mer-m8, (c) 7mer-A1, and (d) 6mer. As shown in Table 3, the seed region types of the miRNAs bta-miR-2418, bta-miR-375, and bta-miR-2468 were all 7mer-m8, except for bta-miR-2479.
Especially, to increase the binding efficiency, we analyzed the G-U wobble base pair parameter using the RNAhybrid software. The binding energies of the miRNAs that bind to the sites with major alleles of the three SNPs were all zero. As a result of the minor alleles shown in Table 3, the binding energy of the miRNA-mRNA seed region with an A allele of SNP g.54853A>G was the lowest at −21.3 kcal mol −1 , whereas the binding energy between the A allele of SNP g.55441A>G and miRNA bta-miR-2479 was the highest at −11.9 kcal mol −1 .

In vivo evidence of the genotype effect of functional SNPs on marbling score
To validate the effect of these functional SNPs on intramuscular fat deposition in vivo, we determined the relationship between three functional SNPs and Korean cattle beef grades (Fig. 2). As shown in Fig. 2, three candidate functional SNPs (g.54853A>G, g.55441A>G, and g.55930C>T) caused a significant difference between beef grades "1 ++ " and "2". As shown in Table 3, the binding affinities of the minor alleles of candidate SNPs were higher than those of their major alleles. Thus, these results suggested that the miRNA-minor alleles of these SNPs act as inhibitors in the regulation of glycerolipid synthesis. As a result of the relationship between these SNPs and beef grades, as shown in Fig. 2, we identified that the "2" beef grade in the group with the minor alleles of these SNPs was significantly lower than that in the group with major alleles (p<0.05).

Discussion
In the present study, we evaluated whether the 3 -UTR of GPAM was influenced by the regulation of glycerolipid synthesis. The 3 -UTR of GPAM is 3689 bp in length and is longer than those of the other genes (Chen et al., 2012). The 3 -UTR has an important influence on the regulation of gene expression because a potential miRNA response element (MRE) is present in this region (Arnold et al., 2012).
Recently, Yu et al. (2017) reported that an SNP in the 3 -UTR of GPAM was significantly associated with fatty acid composition of intramuscular fat and marbling score in a beef cattle population. In the present results, as shown in Tables 1 and 2, we demonstrated that SNPs g.54316A>C, g.54853A>G, g.55441A>G, and g.55930C>T, which are located in the 3 -UTR, were associated with marbling score in the studied Korean cattle steer population. These results coincide with those reported by Yu et al. (2017). Further-32 H. Kim et al.: Determination of the effect of functional SNPs more, as in vivo evidence of the genotype effect of functional SNPs on marbling score, the group with the major alleles of these SNPs had significantly lower binding affinities between the seed region allele and miRNA than those with the minor alleles of these SNP. Thus, these results suggested that these SNPs had an important effect on the regulation of glycerolipid synthesis in the studied Korean cattle steer population. We found that the 3 -UTR of the GPAM gene contained four 15-LOX-DICEs and two CPEs, which regulate mRNA translation and structure stabilization. However, no SNPs were detected in the region containing four 15-LOX-DICEs and two CPEs. In cattle, miRNA targeting has predominantly been associated with the 3 -UTR region of the transcripts derived from ORFs, typically leading to down-regulation through triggering RNA degradation, RNA instability, and/or reduction (He and Hannon, 2004;Li et al., 2011). Thus, as shown in Table 3, four miRNAs were predicted to bind to the seed region including SNPs g.54853A>G, g.55441A>G, and g.55930C>T, using the software programs TargetScan and RNAhybrid. According to the minimum free energy of the binding affinity, each allele of these SNPs had a different binding affinity. Furthermore, this difference in binding affinities led to a significant difference between beef grades "1 ++ " and "2". Thus, our results suggested that the three SNPs, which are present in the 3 -UTR, could be genetic variations influencing the regulation of GPAM gene expression.

Conclusions
Glycerol-3-phosphate acyltransferase, mitochondrial, the enzyme located on the mitochondrial outer membrane, catalyzes the initial and committed step in glycerolipid synthesis in the lipid metabolism of cattle. In the present study, we determined whether the SNPs located in the 3 -UTR of GPAM affect the regulation of gene expression. Of the 11 SNPs detected in the 3 -UTR of GPAM, three SNPs, g.54316A>C, g.54853A>G, g.55441A>G, and g.55930C>T, were significantly associated with marbling score in Korean cattle steer population and showed strong pairwise correlation. Furthermore, we identified four putative miRNAs and found that these SNPs were present in the seed region of these miRNAs. These miRNAs showed a differential binding affinity for each allele of SNPs g.54316A>C, g.54853A>G, g.55441A>G, and g.55930C>T, leading to a significant difference between beef quality grades "1 ++ " and "2". Thus, our results suggested that three SNPs in the 3 -UTR could be genetic variations affecting the regulation of GPAM gene expression.
Data availability. The original data of the paper are available from the corresponding author upon request.
Author contributions. HK and DYO performed the data analyses and wrote the manuscript. YL revised the manuscript and designed the experiment. All authors reviewed and approved the final paper.
Competing interests. The authors declare that they have no conflict of interest.