Genetic Variations and mRNA Expression of Goat DNAH1 and Their Associations with Litter Size

Dynein Axonemal Heavy Chain 1 (DNAH1) encodes proteins which provide structural support for the physiological function and motor structure of spermatozoa (hereafter referred to as sperm) and ova. This study found that three single nucleotide polymorphisms (SNPs), the 27-bp insertion/deletion (InDel) mutations and three exonic copy number variations (CNVs) within DNAH1 were significantly associated with litter size of Shaanbei white cashmere goats (n = 1101). Goats with the wildtypes of these three SNPs had higher litter sizes than other carriers (p < 0.05). II genotype of the 27-bp InDel had the highest litter size compared with ID carriers (p = 0.000022). The gain genotype had the largest litter sizes compared with the loss or medium carriers for the three CNV mutations (p < 0.01). Individuals with the AA-TT-CC-II-M1-M2-M3 and AA-TT-CC-II-G1-G2-M3 combination genotypes had larger litter sizes compared with the other genotypes. This study also showed the DNAH1 expression in mothers of multiple kids was higher than mothers of single kids. These three SNPs, the 27-bp InDel and three CNVs in DNAH1 could be used as molecular markers for the selection of goat reproductive traits.


Introduction
The Shaanbei white cashmere goat (SBWC) excels in the production of both fine cashmere and meat, and therefore plays an important role in the Chinese goat industry. Nevertheless, there is much scope for improvement in their reproductive performance. Poor reproductive performance in goat breeds has impeded the growth of the goat industry [1,2]. Hence, there is a pressing need to develop effective and practical measures for improving reproductive traits in this species [3].
Some animal phenotypes are qualitative traits with high heritability, but most traits associated with reproduction are not. They are often of low heritability, and methods of traditional breeding selection are often not suitable for improving these traits [4]. In contrast, modern molecular labeling techniques can be used to identify beneficial genes, the selection of which will assist in breeding excellent individuals [5]. Some known molecular mutations, single nucleotide polymorphism (SNP) mutations [6], insertion/deletion (InDel) and copy number variant (CNV) variations can be identified quickly and have been employed as useful and efficient molecular markers in the breeding of goats and other livestock. Indeed, SNP, InDel and CNV have been widely used to screen for reproduction-related candidate genes [3,7,8]. Meanwhile, the relationships between SNP, InDel and CNV mutations within certain genes have not been explored yet. Here, we want to explore some mutations within candidate genes and detect whether SNP, InDel and CNV mutations of it are correlated. We also want to identify the mRNA expression of a candidate gene in ewes with single kids and multiple kids. Specifically, we will study the relationship between this gene for reproduction and goat litter size by exploring the associations of SNP, InDel and CNV mutations.
Dynein Axonemal Heavy Chain 1 (DNAH1) encodes proteins which provide structural support for the physiological function and motor structure of spermatozoa (hereafter referred to as sperm) and ova [9,10]. Some genetic mutations within DNAH1 lead to multiple morphological anomalies of the flagella (MMAF) [11], and clinically, MMAF is often accompanied by abnormal sperm morphology and function, leading to male infertility. Many morphological and structural defects can be identified in sperm such as short, curved, irregular morphology and even a complete loss of flagella [12]. Furthermore, abnormal sperm morphology has been documented to affect reproduction and cause infertility in humans and mice [13], dogs [14] and cattle [15]. Meanwhile, the DNAH1 gene is also associated with primary ciliary dyskinesia (PCD) [16,17]. PCD is a type of autosomal recessive inheritance caused by defects in cilia structure or function that are induced by gene mutation. It can seriously affect both male and female reproductive capacity [16]. In males, the sperm tail is a cilia variant so an abnormal structure could cause it to lose motility, which may lead to infertility. In females, the surface of the fallopian tubes is covered by a rhythmically oscillating cilia which ensures the transport of ova; thus, structural dysfunction of cilia can lead to ectopic pregnancy and even infertility [18,19].
As an MMAF-related gene, Cilia and flagella associated protein 43 (CFAP43) gene has been reported to significantly affect kidding traits in goats [20]. Meanwhile, DNAH1 belongs to the MMAF-related gene, and it is a reproductive candidate gene that is associated with both female and male gonadal development and fertility. The DNAH1 gene encodes an inner dynein arm heavy chain that provides structural support between the radial spokes and the outer doublet of the sperm tail, and females can exhibit sub-fertility owing to defective oviduct cilia. Additionally, DNAH1 is also broadly expressed in the testis, ovary and placenta [13,21]. To date, there has been no evidence to confirm whether this gene is related to reproductive traits in goats. Based on the aforementioned studies, we hypothesized that DNAH1 was linked to goat reproductive traits. The objectives of this study were to identify genetic SNP, InDel and CNV polymorphisms in the DNAH1 gene in SBWC goats, as well as explore the difference in DNAH1 expression between mothers of single kids and mothers of multiple kids. We suggest that the results of this study could be used in the design of effective goat breeding programs that would promote growth within the goat industry.

DNA Isolation and Total RNA Extraction form Goat Tissues
For DNA extraction, a total of 1101 adult, female, physically mature SBWC goats were randomly selected from an SBWC farm in Shaanxi Province, China. For RNA extraction, a total of 6 adult female goats were selected, including 3 mothers of single kids and 3 mothers of multiple kids. All does were of similar age and weight and raised under the same external environment [20,22]. After weaning, does were all fed a diet of corn, wheat bran, soybean meal and premix. Farm workers were responsible for the raising and management of all experimental goats and none showed any adverse health conditions. On the SBWC goat farm, veterinarians recorded the number of first litter per goat on a daily basis. The phenotypic distribution of does with different litter sizes is shown in Table 1. Individuals used for sampling were randomly selected from the goat farm to ensure that they were unrelated to each other [7]. According to our previous report, we collected doe ear tissue and ovary samples before the goats were sold or slaughtered. After this sampling step was completed, all samples were placed in an alcohol solution and immediately frozen. The process of sample collection and treatment and the steps of genome-wide extraction were performed as previously reported [20,23]. The quality and purity of each DNA sample was measured by Nanodrop 1000 (Thermo Fisher Scientific Inc., Wilmington, DE, USA). According to previous studies, the mass concentration of DNA extracted from each sample was high [24]. Then, DNA samples were diluted with ddH 2 O to a standard concentration of 20 ng/uL and stored at −20 • C for subsequent genotyping. For RNA assay, in total 6 adult female goats were selected including 3 mothers of single kids and 3 mothers of multiple kids. According to our previous report, we collected doe ear tissue and ovary samples before the goats were slaughtered. Total RNA isolation of ovary and synthesis of cDNA were shown as before [2].

Identification of Candidate Mutations and Primer Design
For SNP detection, using the Ensembl online database (http://asia.ensembl.org/, accessed on 1 December 2020), a pair of primers were designed to amplify fragment containing 10 SNP mutations within the promoter region of DNAH1 gene in 327 goats. These variations including SNP1 (rs665366227), SNP2(rs647320555), SNP3(rs638339874), SNP4 (rs653533484), SNP5 (rs646237158), SNP6 (rs667313924), SNP7 (rs640739292), SNP8 (rs653520370), SNP9 (rs656018411) and SNP10 (rs644545438). For InDel detection, a total of 9 putative InDel mutations of DNAH1 were studied (P1-P9). In accordance with a previous prediction, the corresponding primer sequences were designed according to the 9 mutation positions using Primer software 5.0 (Canada, Premier). For CNV detection, validation of goat exonic CNV1 (g:48593201-48595200), CNV2 (g:48603601-48605200) and CNV3 (g:48617201-48618800) within the DNAH1 gene in chromosome 22 was performed referring to the Animal Omics Database (http://animal.nwsuaf.edu.cn/, accessed on 10 July 2021). Three pairs of primers were designed to amplify these three CNV regions of DNAH1 in 318 goats ( Figure 1). Given that previous studies suggest that there are neither CNVs nor segmental duplication in the goat melanocortin 1 receptor (MC1R) gene, we selected MC1R as the internal reference gene and designed a pair of primers to amplify this target sequence [25]. The specific primer sequences are shown in Table S1.  According to our previous report, we collected doe ear tissue and ovary samples before the goats were sold or slaughtered. After this sampling step was completed, all samples were placed in an alcohol solution and immediately frozen. The process of sample collection and treatment and the steps of genome-wide extraction were performed as previously reported [20,23]. The quality and purity of each DNA sample was measured by Nanodrop 1000 (Thermo Fisher Scientific Inc., Wilmington, DE, USA). According to previous studies, the mass concentration of DNA extracted from each sample was high [24]. Then, DNA samples were diluted with ddH2O to a standard concentration of 20 ng/uL and stored at −20 °C for subsequent genotyping. For RNA assay, in total 6 adult female goats were selected including 3 mothers of single kids and 3 mothers of multiple kids. According to our previous report, we collected doe ear tissue and ovary samples before the goats were slaughtered. Total RNA isolation of ovary and synthesis of cDNA were shown as before [2].

Detection of SNP, InDel and CNV Mutations of DNAH1 Gene
Mutant sites that proved to be polymorphic were further identified in the experimental goats. By the Sanger sequencing, 10 SNP mutations were identified according to whether they were polymorphisms. The sequencing results are shown in Figure 2. Two putative InDel mutations of DNAH1 (DNAH1-P1 and DNAH1-P3) were verified using traditional PCR methods, as previously reported [20,22]. Gel electrophoresis and Sanger sequencing results of the DNAH1-P1 and DNAH1-P3 InDel loci are shown in Figure 3. For the three CNV regions within DNAH1 and the internal reference MC1R, real-time quantitative PCR (qPCR) was used to amplify the target fragment using SYBR Green reactions in triplicate. The reaction volume was 10 µL, and each reaction contained 16 ng of genomic DNA and 5 µL of SYBR Green Mix (BIOER, Hangzhou, China). Thermal cycling conditions were as follows: 95 • C for 10 min followed by 40 cycles of 95 • C for 15 s, 60 • C for 30 s and 72 • C for 20 s.

Detection of SNP, InDel and CNV Mutations of DNAH1 Gene
Mutant sites that proved to be polymorphic were further identified in t mental goats. By the Sanger sequencing, 10 SNP mutations were identified ac whether they were polymorphisms. The sequencing results are shown in Figu putative InDel mutations of DNAH1 (DNAH1-P1 and DNAH1-P3) were verified ditional PCR methods, as previously reported [20,22]

Detection of mRNA Expression of DNAH1 Gene
For the detection of DNAH1 expression assay, the cDNA of both sides of ovary samples including 3 mothers of single kids and 3 mothers of multiple kids were used to explore the mRNA expression of DNAH1 gene by qPCR method. The primers of qPCR and internal gene are shown in Table S1. The procedure was the same as the procedure of CNV method.

Statistical Analyses
As previously reported, the population indexes of the experimental goats were calculated [26], including genotype frequency, allele frequency, homozygosity (Ho), heterozygosity (He), polymorphism information content (PIC), etc. In linkage disequilibrium (LD) analysis, the case of r 2 < 0.33 was considered to represent weak LD; r 2 > 0.33 was considered to represent strong LD; and r 2 = 1 was considered to represent complete LD [27]. The relationships between SNP, InDel and CNV mutations and goat litter sizes were analyzed using chi-square analysis of non-parametric methods and a single-factor analysis method.
In the manuscript, the following fixed model was used: statistical model applied to determine the association of litter sizes with different mutation genotypes:Yijklm = μ + Gi + Bj + Fk + Al+ eijklm, where Yijklm is the phenotypic value of litter size trait, μ is the population mean, Gi is the fixed effect of genotype, Bj represents age effect, Fk is herd effect, Al stands for regional effect and eijklm is the random error [20,28]. Compared with the internal reference gene, the analysis method of copy number variation in each sample: ΔCt = Cttarget gene − Ctreference gene, cycle threshold (Ct) was used to represent the mean of triplicate independent individuals for CNV mutation detection [2,3]. Finally, genotypes of CNV mutations within DNAH1 were divided into three types, including loss (<2), medium (=2) and gain

Detection of mRNA Expression of DNAH1 Gene
For the detection of DNAH1 expression assay, the cDNA of both sides of ovary samples including 3 mothers of single kids and 3 mothers of multiple kids were used to explore the mRNA expression of DNAH1 gene by qPCR method. The primers of qPCR and internal gene are shown in Table S1. The procedure was the same as the procedure of CNV method.

Statistical Analyses
As previously reported, the population indexes of the experimental goats were calculated [26], including genotype frequency, allele frequency, homozygosity (Ho), heterozygosity (He), polymorphism information content (PIC), etc. In linkage disequilibrium (LD) analysis, the case of r 2 < 0.33 was considered to represent weak LD; r 2 > 0.33 was considered to represent strong LD; and r 2 = 1 was considered to represent complete LD [27]. The relationships between SNP, InDel and CNV mutations and goat litter sizes were analyzed using chi-square analysis of non-parametric methods and a single-factor analysis method.
In the manuscript, the following fixed model was used: statistical model applied to determine the association of litter sizes with different mutation genotypes: where Y ijklm is the phenotypic value of litter size trait, µ is the population mean, G i is the fixed effect of genotype, B j represents age effect, F k is herd effect, A l stands for regional effect and e ijklm is the random error [20,28]. Compared with the internal reference gene, the analysis method of copy number variation in each sample: ∆Ct = Ct target gene − Ct reference gene, cycle threshold (Ct) was used to represent the mean of triplicate independent individuals for CNV mutation detection [2,3]. Finally, genotypes of CNV mutations within DNAH1 were divided into three types, including loss (<2), medium (=2) and gain (>2). Experiments were repeated three times and the mean value of intensity ratios ± standard error (SE) was shown [29,30].
The 2 −∆∆Ct method and t-test were used to analyze the expression difference of DNAH1 gene in ovary [31].

Allele and Genotype Frequencies of SNP, InDel and CNV Variations in DNAH1 Gene
The results of allele and genotype frequencies of SNP1 (rs665366227), SNP3 (rs638339874), SNP5 (rs646237158), SNP7 (rs640739292), SNP9 (rs656018411) and SNP10 (rs644545438) are shown in Table 2. The χ 2 test indicated that the genotypic distributions of SNP1, SNP4, SNP5, SNP7, SNP9 and SNP10 were all consistent with Hardy-Weinberg equilibrium (p > 0.05). Agarose gel electrophoresis showed that there were three genotypic variants of the 27-bp InDel (DNAH1-P1): homozygotic insertion, heterozygote and homozygotic deletion type were replaced by II, ID and DD, respectively. In the goat population, we only found a single DD carrier. For the 27-bp InDel site, the "I" allele was associated with the main allele in goat individuals. Meanwhile, the 15-bp InDel (DNAH1-P3) had two genotypes; there were II and ID carriers. Results showed the "D" allele was associated with the dominant allele in goat individuals. The χ 2 test indicated that the genotypic distributions of these two InDels were in Hardy-Weinberg equilibrium (p > 0.05; Table 2). For these three exonic CNVs of DNAH1, the results showed that all three genotypes were present, including loss, medium and gain types.

Association Analysis between SNP, InDel and CNV Variations with Goat Litter Size
Association analysis illustrated that SNP1 (rs665366227), SNP5 (rs646237158) and SNP7 (rs640739292) were significantly associated with goat litter sizes (p < 0.05), and wild genotypes had higher litter sizes than heterozygous genotypes and homozygous mutant genotypes ( Table 3). The 27-bp InDel within DNAH1 was significantly related to goat litter size (p = 0.000022), and the litter sizes of II carriers were higher than those of ID carriers. However, the 15-bp InDel was not related (p = 0.605; Table 3). For exon CNV1, the gain type had significantly larger litter sizes than the medium and loss types. The medium type had significantly higher litter sizes than the loss types (p = 0.000011). For exon CNV2, the gain type had significantly larger litter sizes than the other two types (p = 0.003). For exon CNV3, the gain type also had significantly larger litter sizes than the medium and loss types (p = 0.000187; Table 4).

Relationship between Combination Genotypes and Goat Litter Size
DNAH1-SNP1, DNAH1-SNP5, DNAH1-SNP7, DNAH1-P1, DNAH1-CNV1, DNAH1-CNV2 and DNAH1-CNV3 were all associated with goat litter size (p < 0.05). Next, we aimed to explore the relationship between combinations of genotypes and goat litter size. Among 1101 individuals, we selected combinations of genotypes with sample sizes greater than three for association analysis. There were eight types of combination genotypes (SNP1-  Figure S1). Furthermore, we found that there was strong LD between SNP1, SNP5 and SNP7 mutations (r 2 > 0.33; Figure S2).

mRNA Expression of DNAH1 in Ovary of Goats
The qPCR results proved the mRNA expression of the DNAH1 gene was different in mothers of single kids and mothers of multiple kids. The mRNA expression showed that the DNAH1 gene had higher expression in mothers of multiple kids than mothers of single kids. Although there is no significant difference between them (p > 0.05), the relationship between mRNA expression of goat DNAH1 gene and litter size traits deserve to be further explored. This result indicates that DNAH1 expression might play an important role in litter size traits.

Discussion
The DNAH1 gene is broadly expressed in multiple reproductive tissues, including the testis, ovary and placenta, and is therefore a candidate gene for reproductive traits [21,32]. According to previous reports, mutations within the DNAH1 gene are associated with phenotypes of MMAF and PCD [13,18]. MMAF has been widely studied in male reproductive traits and is often associated with abnormal sperm, which seriously affects male reproductive traits and even leads to infertility [12]. Indeed, as an MMAF-related gene, many papers have identified DNAH1 as playing an essential role in the development of male and female gonads in humans and have explored the clinical applications of DNAH1 [33]. PCD is often accompanied by cilia movement disorder, which is also an autosomal recessive genetic disease. Patients with severe PCD exhibit infertility. Clinically, females can exhibit sub-fertility owing to defective oviduct cilia; on the other hand, males can be infertile because of immotile sperm flagella [17,34,35]. Imtiaz et al. speculated that females may be more sensitive to severe variants of DNAH1 and that these variants may have more severe effects than those identified in males [18]. Specifically, they found that molecular variation in DNAH1 may play a role in PCD, which suggests that its potential contribution to reproduction should be considered in future studies.
Abnormal sperm morphology has been documented to affect reproduction and cause infertility in humans and mice [13], dogs [14] and cattle [15]. We have found that an InDel variation within an MMAF-associated gene, CFAP43 gene, significantly affected litter size in goats [20]. This suggests that MMAF-related genes can also influence reproductive performance in goats. DNAH1 gene encodes an inner dynein arm heavy chain that provides structural support between the radial spokes and the outer doublet of the sperm tail, and females can exhibit sub-fertility owing to defective oviduct cilia. Mice with a homozygous knockout of the orthologous gene are viable but have reduced sperm motility and are infertile.
There are only some reports regarding the relationship between the DNAH1 gene and the reproductive performance of humans and mice. In this study, we found DNAH1 expression of mothers of multiple kids was higher than mothers of single kids (Figure 4). This indicated that the DNAH1 gene might play an essential role in goat productive traits. The current study explored the relationship between several genetic variants and goat litter size to provide scientific guidance for the breeding of goat reproductive traits. In this study, three SNPs, a 27-bp InDel and three exonic CNV loci within DNHA1 were all associated with goat litter size traits (p < 0.05; Table 4). These results suggest that these three SNPs, the 27-bp InDel and the three exonic CNV mutations of DNAH1 would be of great value in the genetic selection of goat reproductive traits. Some studies reported that genetic polymorphisms could affect gene expression and, consequently, some growth and reproduction traits of animals [36,37] and humans [38,39]. For example, there is a single nucleotide mutation in the IGF2 gene of the Bama pig, which can significantly affect skeletal muscle development in this breed. Through genetic breeding of this mutation, the weight traits of this pig breed can be significantly improved [40]. Likewise, we found that these SNP, InDel, and CNV sites within DNAH1 could significantly affect goat reproductive traits, suggesting that these are important variants with potential value in breeding programs owing to their collective action in the modulation of reproductive traits.
Cells 2022, 11, x FOR PEER REVIEW 9 of 12 skeletal muscle development in this breed. Through genetic breeding of this mutation, the weight traits of this pig breed can be significantly improved [40]. Likewise, we found that these SNP, InDel, and CNV sites within DNAH1 could significantly affect goat reproductive traits, suggesting that these are important variants with potential value in breeding programs owing to their collective action in the modulation of reproductive traits. Previous reports identify that linkage relationships between different mutation sites can exert a synergistic effect on the phenotypic traits of livestock [41,42]. Considering that these three SNPs, the 27-bp InDel and the three CNVs simultaneously affect goat litter size, and the relationship between SNP, InDel and CNV mutations has not been previously reported. We speculated on the existence of an association between these SNPs, In-Dels and CNVs in DNAH1. In the current study, we found that these seven genetic mutations all had significant effects on goat litter size. There was a strong linkage disequilibrium between SNP1, SNP5 and SNP7 sites (r 2 > 0.33; Figure S2), and other mutation sites may not affect goat reproduction through strong LD. Thus, the specific mechanisms underlying the functional relationships between these variants deserve to be further studied.
Given that multiple mutations of DNAH1 can play important roles in determining litter size, we believe that these mutations may function synergistically. However, the specific functional mechanism remains to be further explored. Analysis of combined genotypes revealed that individuals with the combination genotypes AA-TT-CC-II-M1-M2-M3 and AA-TT-CC-II-G1-G2-M3 had larger litter sizes relative to other genotypes. This finding illustrates that by selecting for these three SNP sites, the 27-bp InDel and the three CNVs simultaneously, litter size traits could be improved rapidly and efficiently; however, the specific mechanism remains to be further explored.

Conclusions
In this study, three SNPs, the 27-bp InDel and three exonic CNVs within DNAH1 significantly affected litter sizes of Shaanbei white cashmere goats, and DNAH1 expression of ewes with multiple kids was higher than ewes with single kids. Meanwhile, individuals with the combination genotypes AA-TT-CC-II-M1-M2-M3 and AA-TT-CC-II-G1-G2-M3 had the largest litter sizes, which could be used as the effective DNA marker.
Supplementary Materials: The following supporting information can be downloaded at: www.mdpi.com/xxx/s1, Table S1: PCR primers used for detecting in this study; Figure S1: Relationship between the combined genotypes of SNP, the 27-bp InDel, and CNV variations within DNAH1 gene and litter sizes in goats (mean ± standard errors); Figure S2: Linkage disequilibrium analysis of the SNP, InDel and CNV mutations within DNAH1 in Shaanbei white cashmere goats. Previous reports identify that linkage relationships between different mutation sites can exert a synergistic effect on the phenotypic traits of livestock [41,42]. Considering that these three SNPs, the 27-bp InDel and the three CNVs simultaneously affect goat litter size, and the relationship between SNP, InDel and CNV mutations has not been previously reported. We speculated on the existence of an association between these SNPs, InDels and CNVs in DNAH1. In the current study, we found that these seven genetic mutations all had significant effects on goat litter size. There was a strong linkage disequilibrium between SNP1, SNP5 and SNP7 sites (r 2 > 0.33; Figure S2), and other mutation sites may not affect goat reproduction through strong LD. Thus, the specific mechanisms underlying the functional relationships between these variants deserve to be further studied.
Given that multiple mutations of DNAH1 can play important roles in determining litter size, we believe that these mutations may function synergistically. However, the specific functional mechanism remains to be further explored. Analysis of combined genotypes revealed that individuals with the combination genotypes AA-TT-CC-II-M 1 -M 2 -M 3 and AA-TT-CC-II-G 1 -G 2 -M 3 had larger litter sizes relative to other genotypes. This finding illustrates that by selecting for these three SNP sites, the 27-bp InDel and the three CNVs simultaneously, litter size traits could be improved rapidly and efficiently; however, the specific mechanism remains to be further explored.

Conclusions
In this study, three SNPs, the 27-bp InDel and three exonic CNVs within DNAH1 significantly affected litter sizes of Shaanbei white cashmere goats, and DNAH1 expression of ewes with multiple kids was higher than ewes with single kids. Meanwhile, individuals with the combination genotypes AA-TT-CC-II-M 1 -M 2 -M 3 and AA-TT-CC-II-G 1 -G 2 -M 3 had the largest litter sizes, which could be used as the effective DNA marker.
Supplementary Materials: The following supporting information can be downloaded at: https:// www.mdpi.com/article/10.3390/cells11081371/s1, Table S1: PCR primers used for detecting in this study; Figure S1: Relationship between the combined genotypes of SNP, the 27-bp InDel, and CNV variations within DNAH1 gene and litter sizes in goats (mean ± standard errors); Figure S2: Linkage disequilibrium analysis of the SNP, InDel and CNV mutations within DNAH1 in Shaanbei white cashmere goats.