InDel and CNV within the AKAP13 Gene Revealing Strong Associations with Growth Traits in Goat

Simple Summary The AKAP13 gene has been found to be related to bone formation. As an important gene regulating growth and development, whether AKAP13 has any influence on the growth traits of goats is still unknown. Therefore, we hope to explore the InDel and CNV genetic variation of the AKAP13 gene in Shaanbei white cashmere goats and their relationship with growth traits to find effective molecular marker sites. Abstract A-kinase-anchoring protein 13 (AKAP13) is a member of the AKAP protein family that has been found to be associated with bone formation. Thus, we investigated the AKAP13 gene as a potential candidate gene for molecular-marker-assisted selection (MAS) in breeding. Our aim was to explore genetic variations (InDel and CNV) within the AKAP13 gene of Shaanbei white cashmere (SBWC) goats and analyze their relationship with growth traits. Ultimately, we identified three InDel loci (16-bp deletion, 15-bp insertion, and 25-bp deletion) and three CNVs, and the 16-bp and 15-bp loci were significantly associated with goat body length (p < 0.05). Both the 16-bp deletion variant and the 15-bp insertion variant facilitated an increase in body length in goats. In addition to this, there was a certain superposition effect between 16-bp and 15-bp loci, although there was no linkage. Additionally, the CNV1 locus was significantly correlated with body height and body length of goats (p < 0.05), and CNV2 was significantly correlated with chest depth, chest circumference, and cannon circumference of goats (p < 0.05). Individuals with gain type showed excellent growth performance. In conclusion, the InDel and CNV loci that we have identified could possibly serve as effective molecular markers in goat breeding, which is very essential for improving efficiency and success of breeding. Moreover, our findings provide a new avenue for further research into the function of the AKAP13 gene.


Introduction
With the increase in the world population and economic development, the demand for livestock and poultry products has significantly increased [1,2].However, it is difficult to provide an adequate supply of quality goat products due to various factors, such as traditional breeding methods and inefficient production constraints [3][4][5][6].Therefore, in addition to improving feeding methods and management conditions, selecting goat breeds with excellent growth traits is an effective means of promoting the development of the goat industry and improving economic benefits.Shaanbei white cashmere (SBWC) goats are an important and exceptional breed known for their meat and fleece production, and they have been selectively bred in China.SBWC goats are recognized for their rapid growth, resilience to rough feeding conditions, cold tolerance, adaptability, and disease resistance [7].With advancements in molecular biology technology, molecular-markerassisted selection (MAS) breeding has gained momentum and has been applied in industry.Using molecular markers like insertion/deletion (InDel) and copy number variation (CNV) has shown great advantages over traditional breeding [8][9][10].Numerous studies have demonstrated the effective use of InDel and CNV markers in association with growth traits, thereby reducing breeding time and improving outcomes [11].
The rapid progress of high-throughput sequencing technology, along with the continuous improvement of public databases and websites, has made it feasible to identify candidate genes associated with production traits [12][13][14].Herein, A-kinase-anchoring proteins (AKAPs) were chosen because they are a family of proteins with the ability to regulate signal transduction processes [15,16].AKAPs have been found to be widely involved in the regulation of growth, development, metabolism, and reproduction [17][18][19].The main reason for their versatility is that they contain the PKA-binding region and other signaling molecule binding domains, which can bind a variety of signaling molecules, including kinases, phosphatases, adenylate cyclases, GTPases, membrane receptors, and other regulatory proteins, forming a multisignaling complex and regulating their activities [20].AKAP13, as a significant member of this protein family, exhibits similar characteristics.Besides its interaction with PKA, AKAP13 can interact with various protein factors, such as MEK/ERK, p38, and RhoA [21,22].This suggests that AKAP13 is capable of regulating various life processes, including growth.It has been shown that the AKAP13 gene affects heart development in mice, and bone formation and AKAP13 deficiency can lead to osteoporosis and abnormal heart development [23,24].Moreover, the absence of the AKAP13 gene has been linked to weight gain in mice [25].Hence, the AKAP13 gene is a potential candidate associated with growth traits.
While previous studies have verified that AKAP13 may be associated with growth, there is a lack of research on the effect of the AKAP13 gene InDel and CNV variants on growth.Therefore, the objective of this study is to investigate the InDel and CNV loci within the AKAP13 gene and systematically analyze their relationship with growth traits in goats.We aim to identify new molecular markers for MAS breeding in goats and provide foundational data for further understanding the function of the AKAP13 gene.

Materials and Methods
All animal tests performed in this study were conducted under the supervision and guidance of the Faculty of Animal Policy and Welfare Committee of Northwest A&F University, and all procedures were in accordance with the specifications (NWAFU-314020038).

Animal Samples Collection
SBWC goats were kept in similar and suitable conditions at farms in Yulin, Shaanxi, China.A total of 487 ear samples were collected from the adult female goats (2-3 years old) and immediately preserved in liquid nitrogen until used.Various growth traits, including body length, body height, height at hip cross, chest circumference, chest depth, chest width, cannon circumference, and hip width, were measured and recorded for each goat.

Genomic DNA Isolation
Genomic DNA was extracted from the ear tissue samples of the goats using the modified salt-chloroform extraction method [26,27].The concentration and purity of the DNA samples were assessed using Nanodrop 2000 (Thermo Scientific, Waltham, MA, USA).DNA of acceptable quality was used for subsequent experiments.The DNA samples were then diluted to 20 ng/µL with ddH 2 O.In total, 30 samples were randomly selected to construct a DNA mixing pool for preliminary detection of polymorphism at screened InDel sites [28,29].

Primer Design and Genotype Detection
First, we searched the Ensembl database's Variant table (https://asia.ensembl.org/accessed on 4 March 2022) for information on eight InDel loci for the AKAP13 gene in goats.Then, according to the reference sequence of the goat AKAP13 gene (GenBank No: NC_030828.1),eight pairs of primers were designed using NCBI Primer-BLAST (https:// www.ncbi.nlm.nih.gov/tools/primer-blast/accessed on 8 March 2022) and synthesized by Sangon Biotech, China (Table 1).Three CNV loci selected from the Animal Omics Database (http://animal.nwsuaf.edu.cn/accessed on 10 March 2022) were also explored-namely, CNV1, CNV2, and CNV3 with the MC1R gene as the reference gene (Tables 2 and 3).

InDel and CNV Genotyping of AKAP13 Gene
Referring to the previous study, the touch-down PCR procedure was used for amplification.After the end, 6.0 µL PCR products were electrophoresed on a 3.0% agarose gel, which was stained with nucleic acid dyes (Tsingke, Beijing, China) and examined for each individual.Three controls were set for each sample, and then the quantitative real-time polymerase chain reaction (qPCR) results of the AKAP13 gene copy number were analyzed.The qPCR amplification reaction system volume and procedure were as described previously [30], and the result was processed using the method 2 × 2 −∆Ct ; 2 × 2 −∆Ct > 2 was gain type, it < 2 was loss type, and it = 2 was median type.Additionally, the gene structure of the goat AKAP13 gene was drawn as shown in Figure 1.

InDel and CNV Genotyping of AKAP13 Gene
Referring to the previous study, the touch-down PCR procedure was used for amplification.After the end, 6.0 µL PCR products were electrophoresed on a 3.0% agarose gel, which was stained with nucleic acid dyes (Tsingke, Beijing, China) and examined for each individual.
Three controls were set for each sample, and then the quantitative real-time polymerase chain reaction (qPCR) results of the AKAP13 gene copy number were analyzed.The qPCR amplification reaction system volume and procedure were as described previously [30], and the result was processed using the method 2 × 2 −ΔCt ; 2 × 2 −ΔCt > 2 was gain type, it < 2 was loss type, and it = 2 was median type.Additionally, the gene structure of the goat AKAP13 gene was drawn as shown in Figure 1.

Statistical Analyses
The SHEsis program was used to analyze whether these InDel loci were in Hardy-Weinberg equilibrium (HWE).To calculate the genetic parameters, including Ho (observed heterozygosity), He (expected heterozygosity), Ne (effective allele numbers), and PIC (polymorphism information content), the Nei method was used directly [31].Simplified linear models were used to determine the relationship between goat genotypes and each growth trait: where Yij represented the measurement of growth traits for each goat, µ represented the overall mean, Gi represented the effect of genotype, and Eij represented the random error.One-way ANOVA was used to analyze the differences in growth traits between different genotypes and combined genotypes of goats using SPSS software (version 24.0,IBM Corp, Armonk, NY, USA).

Genetic Linkage Analysis
Linkage disequilibrium (LD) analysis was performed on InDel loci using the SHEsis platform [32].D' and r 2 were employed to assess the degree of linkage between InDel loci.

Statistical Analyses
The SHEsis program was used to analyze whether these InDel loci were in Hardy-Weinberg equilibrium (HWE).To calculate the genetic parameters, including Ho (observed heterozygosity), He (expected heterozygosity), Ne (effective allele numbers), and PIC (polymorphism information content), the Nei method was used directly [31].Simplified linear models were used to determine the relationship between goat genotypes and each growth trait: where Y ij represented the measurement of growth traits for each goat, µ represented the overall mean, G i represented the effect of genotype, and E ij represented the random error.One-way ANOVA was used to analyze the differences in growth traits between different genotypes and combined genotypes of goats using SPSS software (version 24.0,IBM Corp, Armonk, NY, USA).

Genetic Linkage Analysis
Linkage disequilibrium (LD) analysis was performed on InDel loci using the SHEsis platform [32].D and r 2 were employed to assess the degree of linkage between InDel loci.

InDel Detection: Genotype Frequency, Linkage Disequilibrium, and Haplotype Analyses of the Goat AKAP13 Gene
All genetic parameters had been calculated and presented in (Table 4).For the P1 InDel locus, the frequencies of the allele "I" and allele "D" were 0.884 and 0.116, respectively.The genotype distribution was not in HWE (p < 0.05).For the P2 InDel locus, DD genotype was much more frequent than II and ID genotypes, and the genotype distribution was in accordance with HWE (p > 0.05).The PIC values indicated that both InDel loci were at low degree of polymorphism (0 < PIC< 0.25).For the P3 InDel locus, it showed medium genetic diversity (0.25 < PIC < 0.5), and "D" allele frequency was higher than "I" allele frequency.The genotype distribution was in accordance with HWE (p > 0.05).

InDel Detection: Genotype Frequency, Linkage Disequilibrium, and Haplotype Analyses of the Goat AKAP13 Gene
All genetic parameters had been calculated and presented in (Table 4).For the P1 InDel locus, the frequencies of the allele "I" and allele "D" were 0.884 and 0.116, respectively.The genotype distribution was not in HWE (p < 0.05).For the P2 InDel locus, DD genotype was much more frequent than II and ID genotypes, and the genotype distribution was in accordance with HWE (p > 0.05).The PIC values indicated that both InDel loci were at low degree of polymorphism (0 < PIC< 0.25).For the P3 InDel locus, it showed medium genetic diversity (0.25 < PIC < 0.5), and "D" allele frequency was higher than "I" allele frequency.The genotype distribution was in accordance with HWE (p > 0.05).The results of the LD analysis are shown in (Figure 6) the D value and r 2 value were calculated (Table 5).These results demonstrated that there was no LD in P1 and P2 InDel loci in this study.Seven haplotypes were detected using haplotype analysis, including haplotype1 (I P1 D P2 I P3 ), haplotype2 (I P1 D P2 D P3 ), haplotype3 (I P1 I P2 I P3 ), haplotype4  6).The haplotype D P1 D P2 I P3 was not found.Following this, we analyzed the association of nine combined genotypes with growth traits in SBWC goats (Table 7).The body length of goats with combined genotype D P1 D P1 -I P2 D P2 was significantly greater than other combined genotypes (p < 0.05).This also demonstrated that these two InDel loci can significantly promote the growth of SBWC goats.
The results of the LD analysis are shown in (Figure 6) the D' value and r 2 value were calculated (Table 5).These results demonstrated that there was no LD in P1 and P2 InDel loci in this study.Seven haplotypes were detected using haplotype analysis, including haplotype1 (IP1DP2IP3), haplotype2 (IP1DP2DP3), haplotype3 (IP1IP2IP3), haplotype4 (IP1IP2DP3), haplotype5 (DP1DP2DP2), haplotype6 (DP1IP2IP2), and haplotype7 (DP1IP2DP2) (Table 6).The haplotype DP1DP2IP3 was not found.Following this, we analyzed the association of nine combined genotypes with growth traits in SBWC goats (Table 7).The body length of goats with combined genotype DP1DP1-IP2DP2 was significantly greater than other combined genotypes (p < 0.05).This also demonstrated that these two InDel loci can significantly promote the growth of SBWC goats.

CNV Detection: Frequency of the Goat AKAP13 Gene Genotypes
According to the 2 × 2 −∆CT method, we divided the CNV types into three classes, including Loss type (0~2), Medium type (2) and Gain type (>2).There are three types (gain, loss, and medium) of CNV1 and CNV2, but only one type of CNV3 (gain) (Table 8).In addition, among the three CNVs, gain type frequency was the highest and was 0.886, 0.692, and 1, respectively.

Association Analysis of Mutations with Growth Traits
The results of the one-way ANOVA showed that only P1 and P2 were significantly associated with body length in SBWC goats (Table 9).For the P1 InDel locus, goats with ID and DD genotypes had significantly higher body lengths than those with II genotype (p = 0.01).This implies that 16-bp deletion had a beneficial effect on goat growth traits.However, for the InDel locus at P2, the body length of goats with the ID genotype was significantly higher than those with the DD genotype (p = 0.014).This means that this 15-bp insertion was a favorable mutation.The superiority of heterozygous genotypes in body length may be due to the complementary effect caused by the comprehensive aggregation of dominant genes in the hybrid offspring.The statistical analyses showed that CNV1 was significantly associated with body height (p = 0.045) and body length (p = 0.046) of goats.For CNV2 locus, there was an association with chest depth (p = 0.011), chest circumference (p = 0.026), and cannon circumference (p = 0.014).For CNV1 and CNV2, the gain type is the higher-frequency type, which had better growth performance than the loss and medium types (Table 10).

Discussion
In this study, we found for the first time that two InDel loci in the AKAP13 gene were significantly associated with growth traits in SBWC goats.Our previous studies have demonstrated that InDel variants in the AKAP13 gene can affect litter size in goats [33].Interestingly, this study also supported our conjecture that the AKAP13 gene could serve as a valid molecular marker for association with growth traits.The three InDel loci we identified belonged to different types, P1 and P3 were deletion mutations, and P2 was an insertion mutation.There was no LD between these three loci.Moreover, based on the PIC values (0 < PIC < 0.25), P1 and P2 were polymorphic to a low degree, indicating low heterozygosity and possibly strong selection [34].The P2 and P3 were in accordance with HWE, while P1 was not.Thus, this suggested that the P1 locus was subjected to a significant degree of environmental or human selection several generations ago, which contributed to the mutations and caused these effects [35].
The results of one-way ANOVA showed that P1 and P2 were significantly associated with goat body length.For the P1 locus, goats with ID and DD genotypes had significantly higher body lengths than the II genotype (p < 0.05), suggesting that the presence of the "D" allele significantly promoted the growth of goats.In addition, we found that this mutation may have a dosage effect.For the P2 locus, the body length of goats with the ID genotype was significantly higher than that of the DD genotype (p < 0.05).While the II genotype was larger than the DD genotype, it did not reach significance (p > 0.05).Therefore, we speculated that the effect of the heterozygous advantage leaded to such a result.The results of the combined genotype analysis showed that goats with the D P1 D P1 -I P2 D P2 genotype exhibited the longest body length.This suggested that although P1 and P2 were not linked, there may be a functional complementarity of these two mutations.
We used the CNVcaller software (https://github.com/JiangYuLab/CNVcaller,accessed on 25 May 2022) to analyze existing whole genome resequencing results and found these three CNVs in the world's goat breeds [36,37], with significant copy number variations in Nubian goat, Longwood goat, and one cashmere goat (Table S1).The resequencing results verified the presence of two loci in Shaanbei white Cashmere goats, after which the population was expanded.The results showed three types (gain, loss, and medium) of CNV1 and CNV2 but only one type of CNV3 (gain).Furthermore, the CNV1 was significantly associated with body height (p = 0.045) and body length (p = 0.046) of goats.For the CNV2 locus, there was an association with chest depth (p = 0.011), chest circumference (p = 0.026), and cannon circumference (p = 0.014).For CNV1 and CNV2, the gain type is the more frequent type, which has better growth performance than the loss and medium types.Thus, this result further validated our previous findings.
AKAP13 is a 10 kDa protein with a PKA-anchoring domain at its N-terminal end, a guanine nucleotide exchange factor (GEF), and a nuclear receptor interaction domain (NRID) at its C-terminal end [38,39].Previous studies have found that AKAP13 is highly expressed in skeletal muscle and bone, suggesting that AKAP13 may play an important role in the regulation of skeletal development [16].Bone formation is a rigorous and complex process throughout animal life, including the regulation of osteoblasts and osteoclasts.Activation of the cAMP/PKA/CREB signaling pathway promotes differentiation of MSCs to osteoblasts, thereby promoting bone growth [40].Thus, AKAP13 can further influence the bone differentiation process by affecting PKA.It is well known that osteoblast differentiation is regulated by various hormones and cytokines, and Rho is one of the most important protein families [25].It is worth noting that recent studies illustrated that a complete RhoA signaling process is required for bone formation [41].AKAP13 can activate RhoA through the GEF structural domain and thus have an important contribution to the differentiation process of osteoblasts [42].Therefore, we hypothesized that variants in the AKAP13 gene may affect the activity of the RhoAGEF structural domain, which in turn affects the activity of RhoA and osteoblast differentiation.On the other hand, runt-related transcription factor 2 (Runx2) is a major regulatory and transcriptional factor for bone development, and the InDel within the Runx2 gene has been found to be associated with reproduction and growth in goats [7,40].Earlier studies have shown that deletion of the AKAP13 gene significantly suppresses its expression [43].Therefore, mutations in the AKAP13 gene may affect bone development by affecting the expression of Runx2.In addition to this, estrogen is also very important for the bone formation process, and lack of estrogen can lead to decreased bone mass and osteoporosis [44].AKAP13 can significantly enhance the ligand-dependent activity of estrogen receptors (ERs) [16].Thus, variation within the AKAP13 gene may affect the binding of AKAP13 to ERs, affecting ER activity and interfering with the function of estrogen.Therefore, mutations in the AKAP13 gene affecting bone growth by affecting the action of estrogen are also one of the possible pathways.
Although both InDel loci we identified were located in the intron region, they can also affect the alteration of animal traits.Mutation in the intronic region might affect the binding of the DNA sequence and DNA binding factors [45,46] as well as the transcriptional efficiency and the stability of mRNA [47,48].Liu et al. (2020) found that the InDel locus located in the intron of the prolactin receptor (PRLR) gene is associated with growth in goats [49].Wang et al. (2020) also found that variants located in the intron region of the IGF2BP1 gene were associated with growth traits in goats [50].Thus, our study further corroborated this conjecture.

Conclusions
In conclusion, in this study, we identified three InDel loci and three CNVs of the AKAP13 gene in SBWC goats and found that the four InDel and CNVs were significantly associated with growth traits and could possibly serve as effective molecular markers for future MAS breeding.

Supplementary Materials:
The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/ani13172746/s1,Table S1: Three CNVs distributions for different goat breeds worldwide.Institutional Review Board Statement: All animal tests performed in this study were conducted under the supervision and guidance of the Animal Welfare Committee of Northwestern Agricultural and Forestry University (NWAFU-314020038-2013-07-10, approved on 10 July 2013), and all procedures were in accordance with their specifications.
Informed Consent Statement: Not applicable.

Figure 1 .
Figure 1.Gene structure of the goat AKAP13 gene.

Figure 1 .
Figure 1.Gene structure of the goat AKAP13 gene.

Funding:
This work was funded by Shaanxi Fundamental Science Research Project for Chemistry and Biology (Grant No. 22JHQ045) and the Shaanxi Province key research and development plan (2022NY-090).

Table 1 .
The primer information of PCR amplification of the AKAP13 gene.

Table 2 .
The primers used for detection of the CNV mutations and the relative expression of the AKAP13 gene.

Table 3 .
Information on the copy number variations within the goat AKAP13 gene.

Table 4 .
Genetic parameters for three InDel loci of the AKAP13 gene.

Table 5 .
D and r 2 value of LD of the AKAP13 gene in SBWC goats.

Table 7 .
The association of goat AKAP13 P1-16 bp and P2-15 bp InDel loci combined genotypes and body growth traits of SBWC goats.

Table 8 .
Typical frequencies of copy number variations within the AKAP13 gene in goat.

Table 9 .
Relationship between three InDel loci of the AKAP13 gene and growth traits in SBWC goats.

Table 10 .
Statistical association analysis of three CNVs of AKAP13 with growth traits in goats.Values with different superscripts within the same line differ significantly at the p < 0.05 (a, b) level.