Identification, Molecular Characterization and Preliminary Functional Analysis of Three SLC27A Full-CDS cDNAs from Water Buffalo (Bubalus bubalis)

Solute carrier family 27 (SLC27A), as known as fatty acid transport proteins (FATPs), are integral transmembrane proteins with long chain fatty acid (LCFA) transport and acyl-CoA synthetase activity (Stahl, 2004). The SLC27As gene family consists of six members, SLC27A1-6, which encode proteins of FATP1-6. SLC27A1, encoding FATP1, was the first member to be identified by screening cDNA that increases LCFA uptake from the cDNA library of 3T3-L1 adipocytes. All FATP members have a highly conserved 311 amino acid tag, called FATP sequence and an AMP binding domain located at the C-terminal. The AMP binding region is responsible for the binding and uptake of LCFA and is commonly found in members of the ACSL family. FATP contains at least one transmembrane domain, the N-terminal is located in the extracellular/intraluminal side and the C-terminal is located in the cytoplasmic side (Anderson and Stahl, 2013; Dourlen et al., 2015). Although SLC27As have sequence similarity, they are differentially expressed in various tissues. SLC27A1 is widely distributed in tissues and organs of Qinchuan cattle and its expression is relatively high in heart, loin muscle and back fat tissues that exhibit rapid fatty acid metabolism, but low in cecum, spleen and lung (Zhao et al., 2015). In chicken, the SLC27A1 transcript had a relatively higher expression in leg muscle, brain and mammary gland (Wang et al., 2011). In humans and mouse, SLC27A2 transcript is expressed primarily in the liver, kidney and lung; SLC27A3 transcript is expressed at lower levels in kidney, liver, heart, lung, ovary and testis; SLC27A4 transcript shows a broader expression pattern with notably high levels in small intestine; SLC27A5 transcript is expressed exclusively in the liver, while FATP6 transcript is expressed primarily in the heart (Mishima et al., 2011). In the mammary gland tissue of Holstein cows during lactation, the mRNA of SLC27A6 was dominant, accounting for 80% of all SLC27As isoforms, followed by the mRNA of SLC27A1, accounting for 13%-50% of all SLC27As isoforms. The abundance of SLC27A5 mRNA was the lowest and SLC27A4 was not detected.
 
So far, the cDNAs of SLC27A1 gene has been isolated and identified in cattle (Ordovás et al., 2006) and this SLC27A1 gene has been mapped in cattle (BTA7q11®q12), sheep (OAR7q11®q12) and goats (CHI5q11®q12) (Ordovas et al., 2005). The sequence of cattle SLC27A2 and SLC27A6 has been annotated in the whole genome database of domestic cattle (Zimin et al., 2009). In addition, studies on SLC27A1, 2 and 6 have been reported in Bovidae species and it has been demonstrated that they play an important role in fatty acid uptake and lipid accumulation. Single nucleotide polymorphisms (g.41534C>T, g.7502A>G, g.5631C>G) in cattle SLC27A1 gene were found to be significantly associated with milk fat content and milk production (Ordovas et al., 2008; Lv et al., 2011). The single nucleotide polymorphisms of g.28470T>C and g.28672G>A were significantly correlated with meat quality traits in Qinchuan cattle (Zhao et al., 2015). Zhang et al., (2021) observed that polymorphism of g.57036072G/T in SLC27A2 gene was closely related to tail lipids. The knockdown of cattle SLC27A6 dramatically downregulated the mRNA abundance of genes associated with FA activation (ACSL4), oxidation (CPT1A) and transport (CD36) in mammary cells and the intracellular triglyceride (TG) content and the ratios of C18:1cis9 and C20:4cis5,8,11,14 were increased, whereas the ratios of C16:0 and C18:0 were decreased (Zhang et al., 2021). It demonstrated that SLC27A6 haplotype was significantly associated with lower saturated fatty acids (SFA) concentrations, lower saturated fatty SFA : unsaturated fatty acids (UFA) ratios and lower lauric acid (12:0) concentrations in bovine milk (Nafikov et al., 2013). At present, few studies have been reported on buffalo SLC27As.
       
Buffalo is an important dairy animal in some parts of the world and the buffalo milk is the second largest source of the world’s dairy industry. Buffalo milk contains more fat, protein, lactose and minerals than cow’s milk and can be used to make butter, good cheeses and other quality dairy products. Buffalo meat is leaner than beef and lower in fat and cholesterol than beef, but the taste is not comparable (Michelizzi et al., 2010; Patel et al., 2018; Singh et al., 2022). Buffaloes are classified into two types: river and swamp buffaloes, with river buffaloes mainly used for milk production and swamp buffaloes mainly used for draught (Mishra et al., 2023). The Binglangjiang buffalo is an indigenous river buffalo found in the Binglangjiang River basin in western Yunnan, China, in 2008. The average milk yield of this buffalo during lactation was 2452.2±553.8 kg, with protein and milk fat percentage of 4.60% and 6.82%, respectively. To date, the genes associated with milk production traits in buffaloes and their mechanisms of action are unknown. In this study, the full-length CDSs of SLC27A1, 2 and 6 were isolated from Binglangjiang buffalo and bioinformatic analyses of the physicochemical characteristics, structure and function of the protein encoded by these three genes. The results here will lay a foundation for further understanding the biological function of these three genes in water buffalo.

Sample collection, RNA extraction and cDNA synthesis
 
The experiment was conducted from December 2023 to May 2024 at the College of Chemistry, Biology and Environment, Yuxi Normal University. Three mammary tissues of Binglangjiang buffalo were obtained from Binglangjiang buffalo breeding farm in Tengchong city, Yunnan Province, China. Total RNA of mammary tissue was extracted by RNAiso Plus (TaKaRa, Dalian, China). The concentration and purity of RNA were detected by NanoDrop 2000 UV-Vis spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA) and the integrity of total RNA was determined by 1.5% agarose gel electrophoresis. The PrimeScript TM 1^st Strand cDNA Synthesis Kit (TaKaRa, Dalian, China) was used to synthesize cDNA.
 
Isolation of buffalo SLC27A1, SLC27A2 and SLC27A6 genes
 
Three specific primer pairs were designed based on the predicted mRNA sequence of buffalo SLC27A1 (accession no. XM_006057903.3), SLC27A2 (accession no. XM_006077795.4) and SLC27A6 (accession no.XM_0060 77795.4) using Primer Premier 5.0, respectively (Table 1). The complete coding region (CDS) of buffalo SLC27A1, SLC27A2 and SLC27A6 gene was amplified using the cDNA of mammary tissues. PCR system was 50 μL containing 1 µL Pfu DNA Polymerase (2.5 U/µL, TIANGEN BIOTECH (Beijing) CO., LTD), 1 μL (10 μmol/L) each primer, 1 μL 10×Pfu Buffer, 4 μL dNTP Mixture (2.5 mM), 2 μL cDNA and 40 μL sterile water. PCR was executed with an initial step at 94°C for 5 min, then 35 cycles of 94°C for 30 s, annealing for 30 s and 72°C for 3 min and finally 72°C for 5 min. The PCR products were performed by 1.5% agarose gel electrophoresis and were subsequently purified with Gel Extraction Kit (OMEGA). The purified PCR products were ligated into pEASY®-T1 Simple Cloning Kit (TransGen Biotech, Beijing, China) and cloned into Trans1-T1 Phage Resistant Chemically Competent Cells (TransGen Biotech, Beijing, China). Then the target DNA fragment of recombinant plasmids were sequenced bidirectionally by Shanghai Biological Engineering Technology Services Co., Ltd. (Shanghai, China).
 

 
Bioinformatics analysis
 
The open reading frame (ORF) of the cloned sequences for SLC27A1, 2 and 6 genes were analyzed using ORF Finder (http://www.ncbi.nlm.nih.gov/gorf/). Then, based on the sequence of ORF as the query sequence, the homologous search was performed using the BLAST program (https://blast.ncbi.nlm.nih.gov/Blast.cgi). The identity and divergence of amino acid sequences were built by MegaAlign program of Lasergene 7 software package (DNAStar, Inc., USA) and further visualized by TBtools software. The physicochemical characteristics, structure, subcellular localization, interaction protein and functional analysis of buffalo FATP1, 2 and 6 proteins were obtained by online softwares that are shown in Table 2.
 

       
Genome GTF files of SLC27A1, 2 and 6 genes for different species of Bovide were downloaded from NCBI Datasets (https://www.ncbi.nlm.nih.gov/datasets/). The mRNA and UTR information were obtained by TBtools software and then the gene transcriptional structure was visualized in Gene Structure Display Server 2.0 (http://gsds.gao-lab.org/). A neighbor-joining (N-J) phylogenetic tree was built by Mega 7 according to amino acid sequences using bootstrap method and Jones-Taylor-Thornton (JTT) model (Kumar et al., 2016). The amino acid sequences of FATP1, 2 and 6 for each species were submitted to the MEME website (http://meme-suite.org/tools/meme) and the NCBI Batch Web CD-Search tool (https://www.ncbi.nlm.nih.gov/Structure/wrpsb/bwrpsb.cgi) to obtain their motif structures and conserved domains information, respectively. Finally, the TB-tools software was used to integrate and visualize the analysis results of motifs and conserved domains.

Isolation and identification of buffalo SLC27A1, 2 and 6 genes
 
The ORF Finder program analysis showed that the sequences obtained from buffalo SLC27A1, 2 and 6 contained 1941 bp, 1863 bp and 1905 bp ORFs, respectively. And the homology search showed that the ORFs had high identity with that SLC27A1, 2 and 6 of other species of Bovidae (more than 97%). Therefore, the ORFs were identified as the CDSs of buffalo SLC27A1, 2 and 6 . The CDS base composition of A, G, T and C for SLC27A1 gene was 17%, 33.49%, 18.19% and 31.32%, respectively and the content of G+C was 64.81%. The CDS base composition of A, G, T and C for SLC27A2 was 26.63%, 25.44%, 25.39% and 22.54%, respectively and the content of G+C was 47.99%. The CDS base composition of A, G, T and C for SLC27A6 was 28.66%, 24.36%, 27.87% and 19.11%, respectively and the content of G+C was 43.46%. The CDSs and their deduced amino acids are presented in Fig 1.
 

       
The transcriptional structure of buffalo SLC27A1, 2 and 6 and their comparison with those of other Bovidae species are e showed in Fig 2. There was only one transcript of SLC27A1 in buffalo (Bubalus bubalis), yak (Bos mutus) and bison (Bison bison bison) with 13 exons in buffalo and yak while 15 exons in bison. However, the SLC27A1 of cattle (Bos taurus), goat (Capra hircus) and sheep (Ovis aries) contain 3, 4 and 4 transcripts, respectively and they were all consist of 13 exons. The main differences between these transcripts of SLC27A1 in these bovine species were at the 5' UTR end and 3' and their CDSs were all the same length of 1941 bp, encoding 647 amino acid residues. Especially, among the six transcripts from zebu (Bos indicus), transcripts XM_019963878.1 and XM_019963877.1 had one more exon in front of the first exon than the others, but lacking exon 2, which resulted in their CDSs being 1,557 bp and 1,563 bp in length and encoding 519 and 521 amino acid residues, correspondingly. For SLC27A2 gene, buffalo, cattle, zebu and sheep had only one transcript, but bison and yak contained two transcripts. The transcript of XM_005892159.1 lacked exon 3 (CDS length: 1704 bp, encoding 568 amino acids), others for all Bovidae species contain 10 exons, with a CDS length of 1863 bp, encoding 620 amino acid residues. The SLC27A6 gene of buffalo, yak, goat and sheep contained only one transcript, while cattle and zebu contained two transcripts. Interestingly, compared to buffalo, goat and sheep, the bovine genera (cattle, zebu and yak) had a 3-base insertion before the 28^th base of the CDS, resulting in their coding product being one amino acid longer than that of the other species (the length of CDS: 1908 bp, amino acids: 635 aa). In particular, the transcript XM_005209093.4 of cattle lacked exon 5 than others and had a CDS length of 1713 bp, encoding 571 amino acids.
 

 
Sequence similarity and phylogenetic analysis
 
The alignment of amino acid sequences showed that buffalo FATP1, 2 and 6 have high consistency with those of other species, especially the species of Bovidae (Fig 3). The sequences of buffalo FATP1, 2 and 6 all have more than 88.6% identity with that of other species of Bovidae and more than 76.4% identity with other mammalian species (Fig 3). The phylogenetic tree based on the AASs showed that FATP1, 2 and 6 were clustered in their own clade and the genetic relationship between buffalo and other species of Bovidae was relatively close (Fig 4).
 

 

       
Motifs and conserved domains were analyzed based on the amino acid sequences of the above mentioned bovine families (Fig 5). In total, of the 10 motifs identified, all were present in the FATP1, A2 and A6 proteins of the bovine species with the exception of motifs 6, 8 and 10. The FATP1 of zebu lacked motifs 6 and 10. It’s worth noting that in all bovine species, FATP2 and A6 lacked the motif 8. All of the FATP1 proteins contained a AFD_class_I superfamily, while FATP2 and A6 contain hsFATP2a_ACSVL_like.
 

       
The basic characteristics of FATP1, 2 and 6 proteins of buffalo and cattle (FATP1: XP_005208560.1; FATP2: NP_001179792.1; FATP6: NP_001094639.1) were compared in Table 3. The results showed that the basic molecular properties of FATP1, 2 and 6 were similar between buffalo and cattle. The values of GRAVY for FATP1, 2 and 6 indicated FATP1 was hydrophobic while FATP2 and FATP6 were hydrophilic proteins. The instability index (II) for buffalo FATP1, 2 and 6 were less than 40, indicating that they were all stable proteins. The FATP1, 2 and 6 contained a transmembrane helix, but had no signal peptide for both buffalo and cattle.
 

       
The secondary structure of buffalo FATP1 consist of 32.97% α-helix (213 AAs), 23.37% extended strand (151 AAs), 8.67% β-turn (56 AAs) and 34.98% random coil (226 AAs); Buffalo FATP2 contained 38.23% α-helix (237 AAs), 20.32% extended strand (126 AAs), 8.55% β-turn (53 AAs) and 32.90% random coil (204 AAs) and buffalo FATP6 was composed of 38.64% α-helix (245 AAs), 19.72% extended strand (125 AAs), 8.20% β-turn (52 AAs) and 33.44% random coil (212 AAs) (Fig 6). The predicted 3D structure displayed that the sequence identity between buffalo FATP1 and human (template: Q6PCB7.1. A) was 92.41% and the coverage rate were 100%. Buffalo FATP2 showed 85.48% sequence identity with human and 100% coverage (template: O14975.1. A). For FATP6, the sequence consistency between buffalo and Rattus norvegicus was 76.82% and the coverage was 97% (template: D4A2B8.1. A) (Fig 7).
 

 

 
Signal peptide and transmembrane helix
 
Prediction showed that buffalo FATP1, 2 and 6 had no N-terminal signal peptide, which indicated that these proteins were non-secretory protein. Transmembrane analysis by DeepTMHM demonstrated that FATP1, 2 and 6 all contained a transmembrane domain (amino acids of 14-41, 4-19 and 17-34, respectively) (Fig 8).
 

 
Subcellular location and molecular function
 
Predictions of subcellular localization indicated that buffalo FATP1, 2 and 6 were functionally in the inner mitochondrial membrane, endoplasmic reticulum membrane and cell membrane, respectively. Buffalo FATP1, 2 and 6 were engaged in many biological processes such as adiponectin-activated signaling pathway (GO:0033211), phosphatidylethanolamine biosynthetic process (GO:0006646), long-chain fatty acid import into cell (GO:0044539), phosphatidylinositol biosynthetic process (GO:0006661), medium-chain fatty acid transport (GO:0001579), response to insulin (GO:0032868), phosphatidylcholine biosynthetic process (GO:0006656), phosphatidylserine biosynthetic process (GO:0006659), negative regulation of phospholipid biosynthetic process (GO:0071072), phosphatidic acid biosynthetic process (GO:0006654) and cardiolipin biosynthetic process (GO:0032049). Their molecular functions were mainly fatty acid transmembrane transporter activity (GO:0015245), long-chain fatty acid transporter activity (GO:0005324), long-chain fatty acid-CoA ligase activity (GO:0004467).
 
Amino acid modifications
 
In this study, five, seven and five putative types of modification sites were found in buffalo FATP1, 2 and 6 by Prosite Scan, including N-myristoylation sites, protein kinase C phosphorylation sites, tyrosine kinase phosphorylation site 1, cAMP-and cGMP-dependent protein kinase phosphorylation sites, N-glycosylation sites, casein kinase II phosphorylation sites and amidation sites (Table 4). The buffalo FATP1 and FATP6 did not contained tyrosine kinase phosphorylation sites and amidation sites.
 

       
The CDSs of buffalo SLC27A1, 2 and 6 were cloned for the first time, with lengths of 1941 bp, 1863 bp and 1905 bp, respectively. Their CDSs exhibited more than 97% identity to those of bovine species such as cattle and goat and their amino acid sequences demonstrated more than 88.6% identity. Phylogenetic tree analyses demonstrated that buffalo SLC27A1, 2 and 6 were clustered with bovine species, indicating that these three proteins are functionally more similar to bovine species. Studies on cattle revealed that polymorphisms in the SLC27A1 gene sequence were significantly correlated with meat quality and milk fat content (Lv et al., 2011; Zhao et al., 2015). Studies on the FATP1 protein showed two opposite function for this gene: inhibitory and promotional effects on lipid deposition. Huang et al., (2022) demonstrated that FATP1 promotes the adipogenic differentiation of buffalo muscle derived cells, but did not affect or even inhibit adipogenic differentiation of adipose-derived cells. Some studies showed that FATP1 enhanced lipid accumulation in cells derived from muscle such as quail QM-7 muscle cells, human muscle cells and porcine intramuscular preadipocytes (Garcia-Martinez et al., 2005; Chen et al., 2017; Qiu et al., 2017), while others have reported that FATP1 does not affect, or even inhibits lipid accumulation in muscle tissue. Muscle contraction requires continuous energy expenditure and FATP1 promotes fatty acid oxidation to provide energy. In rat skeletal muscle and cardiomyocytes, FATP1 was predominantly localized in mitochondria, thereby mediating fatty acid oxidation in these tissues. Mitochondrial localization of FATP1 in muscle tissue correlates with fatty acid oxidation by FATP1 (Nickerson et al., 2009; Sebastian et al., 2009; Huang et al., 2021). Accordingly, it is hypothesized that the opposites of FATP1 functions may be linked to its distinct tissue and subcellular localization. The bioinformatics analyses conducted in the present study indicate that the buffalo FATP1, 2 and 6 were localized in inner mitochondrial membrane, endoplasmic reticulum membrane and cell membrane, respectively. Further investigation is required to determine their specific functions in different tissues and cells of buffalo.
       
Previous studies have found FATPs contained one or more membrane-spanning regions, multiple membrane-associated regions and an AMP-binding motif in the intracellular region. And fatty acids transport through FATP was ATP-dependent, the AMP-binding region was required for transport. Loss-of-function mutations in the AMP-binding domain prevented FA uptake, suggesting AMP-binding was directly involved in or coupled to the transport mechanism (Lewis et al., 2001; Wade et al., 2021). In this study, the conservative domain analyses revealed that buffalo FATP1 contained AMP-binding domain and belonged to AFD_class_I superfamily, FATP2 and A6 contained the domain of AMP-binding and hsFATP2a_ACSVL_like and belonged to AFD_class_I superfamily. The AFD_class_I superfamily includes acyl- and aryl-CoA ligases, as well as the adenylation domain of no ribosomal peptide synthetases and firefly luciferases. The adenylate-forming enzymes catalyze an ATP-dependent two-step reaction to first activate a carboxylate substrate as an adenylate and then transfer the carboxylate to the pantetheine group of either coenzyme A or an acyl-carrier protein. The above findings indicating that the buffalo FATP1, A2 and A6 protein functions as a fatty acid transporter through the AMP-binding domain. Furthermore, this study revealed that the N-terminus of the buffalo FATP1, FATP2 and FATP6 proteins all contained a transmembrane helix, which was consistent with studies in species such as human and mouse (Anderson and Stahl, 2013).
       
After downloading the GTF files of each bovine species from the NCBI database and analyzing the structure of their transcription regions, we found that buffalo SLC27A1, 2 and 6 have only one transcript and the cloning result also showed that their CDSs were not found selectively splicing. However, there were selective splicing in other bovine species including 5' exon hopping, variable 5' splicing site, variable 3' splicing site. Previous studies have shown that alternative splicing can have a variety of functional effects on its encoded products (Florea, 2006). In addition, two isoforms of FATP2 have been identified in humans: FATP2a and FATP2b, of which FATP2b lacks acyl-coenzyme synthase activity due to the deletion of exon 3, but is fully functional in fatty acid transport (Melton et al., 2011). These selective splicing effects on the function of FATP1, 2 and 6 in each species deserve further in-depth study. Several modification sites were found in buffalo FATP1, 2 and 6, including N-myristoylation sites, protein kinase C phosphorylation sites, tyrosine kinase phosphorylation site 1, cAMP-and cGMP-dependent protein kinase phosphorylation sites, N-glycosylation sites, casein kinase II phosphorylation sites and amidation sites. Whether these putative modification sites play crucial roles in the process of fatty acid transport needs further investigation.

In this study, the SLC27A1, 2 and 6 gene were isolated and characterized in Binglangjiang buffalo. The results indicate that the molecular characteristics and functions of buffalo SLC27A1, 2 and 6 are very similar to those of other mammals, especially species of Bovidae. The proteins of buffalo FATP1, 2 and 6 do not have a signal peptide but have a transmembrane structure and were localized to various organelle membranes and cell membranes, which may functionally in the transport of fatty acids through AMP binding sites.

Fangting Zhou was responsible for obtaining the experimental data, writing, revising and approving this manuscript.

The original data of the paper are available from the corresponding author upon request.

This study did not damage the welfare of animals. All efforts were made to minimize any discomfort during the experiment. All sample collection procedures adhered to the “Laboratory animals-General requirements for animal experiment (GB/T-35823-2018) and directive 2010/63/EU on the protection of animals used for scientific purposes.

This research has been supported by the Special Basic Cooperative Research Programs of Yunnan Provincial Undergraduate Universities, China (grant no. 202101BA070001-072) and Yunnan Fundamental Research Projects, China (grant no. 202201AU070175).

The authors declare that they have no conflict of interest.