Next Article in Journal
Management of Donkeys in Assisted Interventions: A Snapshot
Previous Article in Journal
Effects of Probiotic Lactiplantibacillus plantarum HJLP-1 on Growth Performance, Selected Antioxidant Capacity, Immune Function Indices in the Serum, and Cecal Microbiota in Broiler Chicken
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Genome-Wide Analysis of Milk Production Traits and Selection Signatures in Serbian Holstein-Friesian Cattle

1
Department of Biology, Faculty of Veterinary Medicine, University of Belgrade, Bul. Oslobodjenja 18, 11000 Belgrade, Serbia
2
Department of Animal Science, Biotechnical Faculty, University of Ljubljana, Groblje 3, 1000 Ljubljana, Slovenia
3
Department of Reproduction, Fertility and Artificial Insemination, Faculty of Veterinary Medicine, University of Belgrade, Bul. Oslobodjenja 18, 11000 Belgrade, Serbia
*
Author to whom correspondence should be addressed.
Animals 2024, 14(5), 669; https://doi.org/10.3390/ani14050669
Submission received: 5 January 2024 / Revised: 5 February 2024 / Accepted: 14 February 2024 / Published: 21 February 2024
(This article belongs to the Section Animal Genetics and Genomics)

Abstract

:

Simple Summary

This study focuses on identifying the genetic factors associated with milk-related traits in dairy cows. This research employed various genomic techniques to discover a significant association between a specific gene marker and milk protein concentration. The study also assessed inbreeding levels and revealed insights into the genetic mapping of dairy cows. We identified a marker that is significantly associated with milk protein concentration in first lactation (adjusted to 305 days) and, in addition to this marker, we also revealed genomic regions under selection pressure for other economically important traits. Moreover, we revealed low inbreeding levels among the tested animals. These findings contribute to enhancing breeding programs for Holstein-Friesian cattle and improving milk production.

Abstract

To improve the genomic evaluation of milk-related traits in Holstein-Friesian (HF) cattle it is essential to identify the associated candidate genes. Novel SNP-based analyses, such as the genetic mapping of inherited diseases, GWAS, and genomic selection, have led to a new era of research. The aim of this study was to analyze the association of each individual SNP in Serbian HF cattle with milk production traits and inbreeding levels. The SNP 60 K chip Axiom Bovine BovMDv3 was deployed for the genotyping of 334 HF cows. The obtained genomic results, together with the collected phenotypic data, were used for a GWAS. Moreover, the identification of ROH segments was performed and served for inbreeding coefficient evaluation and ROH island detection. Using a GWAS, a polymorphism, rs110619097 (located in the intron of the CTNNA3 gene), was detected to be significantly (p < 0.01) associated with the milk protein concentration in the first lactation (adjusted to 305 days). The average genomic inbreeding value (FROH) was 0.079. ROH islands were discovered in proximity to genes associated with milk production traits and genomic regions under selection pressure for other economically important traits of dairy cattle. The findings of this pilot study provide useful information for a better understanding of the genetic architecture of milk production traits in Serbian HF dairy cows and can be used to improve lactation performances in Serbian HF cattle breeding programs.

1. Introduction

The Holstein-Friesian (HF) is the leading dairy breed in the world; its genetic potential for high milk yields is being exploited in virtually all developed countries. Milk production and milk composition are the most carefully selected traits in dairy cattle breeding programs [1,2]. The discovery of molecular markers associated with the aforementioned traits and the application of marker-assisted selection (MAS) made a huge impact on dairy cattle breeding in recent decades. To improve the genomic evaluation of milk-related traits and to comprehend the underlined molecular mechanisms, it is essential to identify the trait-associated genomic areas and candidate genes [3]. When introduced, genome-wide association studies (GWASs), took the lead over traditional QTL methods both in power of detecting causative gene variants and in defining narrower genomic regions containing causal variants [4]. Since GWASs are ideal for discovering genes that encode complex traits and their mechanisms [5], they have been widely conducted in recent years to identify relationships between genomic variants and economically significant traits in dairy cattle populations. As a result, a large number of candidate genes and quantitative genomic areas were detected to be associated with economically important traits, including milk-related traits, fertility, growth, and meat and carcass quality [3,4,6,7,8,9,10,11,12,13,14,15,16,17]. GWASs conducted in various cattle breeds utilized different SNP chip platforms and were focused on the identification of potential candidate genes and genomic regions responsible for milk production traits including milk yield (MY), fat yield (FY), and protein yield (PY). The most frequently reported candidate genes (DGAT1, GHR, MAPK15) and SNPs associated with milk production traits are located on chromosome 14 [4]. Other genes of interest harbored on chromosomes 14 and 20 are GNA14, PTBP2, and U6 which were closely associated with a higher MY [6,12,18], respectively. ITPR2 [19], ABCC9 [20,21], CPSF1 [21], PDE4 [13], and METTL15 [19] genes were reported to influence the expression of FY, while the genes associated with PY were PDHA2 [19], CTBP2 [19], MAPK9 [22], HPS3 [23], ARFGEF [24], SLCO1A2 [21], and MFSD1 [13].
The availability of a large number of known single nucleotide polymorphisms (SNP) contributed to the research regarding genomic selection and genomic inbreeding [25,26,27,28]. Genomic selection based on SNP chip genotypes relies on linkage disequilibrium (LD) between the QTL and the SNPs [29]. Traditionally, inbreeding was calculated by estimating inbreeding coefficients based on the pedigrees. The progress of modern genotyping methods enabled the increased use of genomic information for more precise estimation of inbreeding coefficients [30]. Several studies have shown that the characterization of inbreeding based on long runs of homozygous (ROH) genotypes presents a more accurate method of measuring individual autozygosity compared to the estimation of total inbreeding based on pedigree data [31,32,33]. The accuracy of the pedigree-based estimation of the inbreeding coefficient depends on the quality and completeness of the pedigrees. Pedigrees provide information on the expected proportion of the genome that is identical by descent (IBD), while the quantification of realized IBD can be measured directly with molecular markers. Furthermore, the genomic inbreeding coefficient provides an assessment of inbreeding by measuring the actual homozygosity. Moreover, the analysis of runs of homozygosity (ROH) provides insights into historical inbreeding patterns within the population and reveals genomic regions with inbreeding effects. Analyses of ROH regions can be used for animal breeding plans in order to minimize the adverse effects of inbreeding. Several studies have investigated the association of ROH with unfavorable gene variants in farm animals [34,35,36]. The signatures of natural and artificial selection “imprinted” on the genome can be traced back and contribute to a better understanding of the evolutionary processes that shaped the cattle genome [37,38,39]. The identification of genes and genomic regions affected by selection is essential to understand the biological mechanisms underlying the phenotypic differences observed between different species of farm animals that have different purposes and are influenced by different environmental conditions [37].
To the best of our knowledge, there is no official HF breeding program in Serbia, and in past decades, we were dependent on “genetics from abroad”, which is one of the reasons why we conducted this study. Furthermore, SNP chip analyses have never been performed before on Serbian HF dairy cattle, which is why this presents a very significant pilot study. Therefore, we decided to perform a GWAS using different software to analyze the association between SNPs and important production traits in a small population of Serbian HF cattle. In addition, the association between SNPs and inbreeding was assessed.

2. Material and Methods

2.1. Sampling

The material used for DNA extraction and further molecular-genetic analyses was collected from 334 chosen cows of the HF breed (that just finished their second lactation) from the farm units of the Al Dahra Corporation, Belgrade, Serbia (44.928635, 20.438109), the largest farm in Serbia, with approximately 4000 lactating cows. All animals were in good condition, without significant health problems within the first two lactations, and all were kept under the same conditions, in a tied system and received the same feed. All of these were criteria for inclusion in the study.
Blood was collected aseptically from the coccygeal vein in 10-mL vacutainers with added anticoagulant (EDTA) and stored at 4 °C until DNA extraction. All blood samples were collected during regular activities related to the implementation of animal health measures. New sterile gloves, needles, and vacutainers were used for each sample collection. All cows were kept under the same conditions, in a tied farm system and fed the same diet.
A total of 9289 AT4 milk records were taken from the database of the farm units of Al Dahra Corporation. Phenotypic data include total milk yield, milk protein, milk fat concentration, as well as milk yield, milk protein and milk fat concentration adjusted to 305 days for 334 genotyped HF cows. Cows were in their first to second parity. Analyzed traits were: trait 1: total milk yield in lactation 1, trait 2: total milk yield in lactation 2, trait 3: total milk protein concentration in lactation 1, trait 4: total milk protein concentration in lactation 2, trait 5: total milk fat concentration in lactation 1, trait 6: total milk fat concentration in lactation 2, trait 7: milk yield adjusted to 305 days in lactation 1, trait 8: milk yield adjusted to 305 days in lactation 2, trait 9: milk protein concentration adjusted to 305 days in lactation 1, trait 10: milk protein concentration adjusted to 305 days in lactation 2, trait 11: milk fat concentration adjusted to 305 days in lactation 1, trait 12: milk fat concentration adjusted to 305 days in lactation 2.

2.2. DNA Extraction and Genotyping

DNA extraction was performed using a commercial kit for BLOOD Version II (Lucigen Corporation, St. Middleton, WI, USA) according to the protocol of the Animal Production and Health Laboratory, International Atomic Energy Agency—IAEA, Seibersdorf, Vienna, Austria, which was additionally adapted for the purposes of this study. Briefly, the extraction was conducted from 1 mL of whole blood and consisted of removal of red blood cells, lysis of white blood cells, and precipitation and pelleting of the genomic DNA. The extracted DNA samples were quantified using a UV-VIS spectrophotometer (BioSpec-nano, Shimadzu Scientific Instruments, Kyoto, Japan) in order to check the DNA quality.
After DNA extraction, we performed genotyping of the extracted samples using the SNP chip Axiom Bovine BovMDv3 (ThermoFisher Scientific, Waltham, MA, USA) with 63,648 markers in collaboration with the IAEA and their laboratory (Animal Production and Health Laboratory, IAEA, Seibersdorf, Vienna, Austria). The analyses were performed according to the established protocol of the above-mentioned laboratory.

2.3. Quality Control of SNP Chip Data

Quality control of SNP chip data was performed using SNP & Variation Suite v8.9.1 software (Golden Helix, Inc., Bozeman, MT, USA). Samples with a genotyping success rate ≤ 0.90 were removed. In addition, markers were removed if the genotyping success rate was <0.95, if they had >2 alleles, or if the minor allele frequency (MAF) was <0.01. The data were then filtered to remove markers that were in linkage disequilibrium (LD), while leaving certain markers to represent groups that were linked to each other. This procedure left 29,503 out of 63,648 markers in 307 animals. MAF and LD purification was skipped to analyze the ROH segments, leaving 52,934 markers in 307 animals. In the data set for the ROH analyses (52,934 markers), we did not include markers from sex chromosomes, while in the GWAS analysis (29,503 markers), we included markers from the X chromosome in addition to the markers from autosomes.

2.4. GWAS

A GWAS was performed using the GAPIT3 library for R [40] with GLM (general linear model), MLM (multiple loci mixed model), FarmCPU (fixed and random model circulating probability unification), Super (settlement of MLM under progressively exclusive relationship), ECMLM (enriched CMLM), CMLM (compressed MLM), and Blink (Bayesian information and linkage-disequilibrium iteratively nested keyway). We considered the population structure and genetic relations within the studied population by performing principal component analysis and calculating the genotype relatedness matrix (kinship). The kinship matrix and the first three principal components (PCs) were fitted as covariate variables in the association models. We used the Bonferroni adjustment to determine GWAS significance cut-off (α  =  0.01).

2.5. Identification of ROH Segments

We identified ROH segments using SNP & Variation Suite v8.9.1 software (Golden Helix, Inc., Bozeman, MT, USA). We defined ROH as 25 or more consecutive homozygous markers on a segment with a minimum length of 1000 kilo base pairs (Kbp). Heterozygous markers were not allowed and there were no more than five unread SNPs. ROHs were classified based on their length (1–2 Mbp, 2–4 Mbp, 4–8 Mbp, 8–16 Mbp and >16 Mbp).

2.6. Genomic Inbreeding Coefficient Based on ROH

The genomic inbreeding coefficient (FROH) [41] was calculated for each animal as the proportion of autosomal genome in ROH. The autosomal genome of Bos taurus is set at 2,512,082,506 bp. We calculated (FROH) including all classes of ROH segments (1–2 Mbp, 2–4 Mbp, 4–8 Mbp, 8–16 Mbp and >16 Mbp). By adding the ROH segments of each class, we also calculated FROH1–2 Mb, FROH2–4 Mb, FROH4–8 Mb, FROH8–16 Mb and FROH > 16 Mb.

2.7. Detection of ROH Islands

We estimated the frequency of SNPs in ROH (%) and then plotted them against their positions along the chromosomes (Manhattan plot) using the statistical software library for R [42]. The minimum threshold for ROH island detection was set at 30%, meaning that ROH must be present in at least 30% of the population for a locus to be included in an ROH island [43].

3. Results

3.1. Genome-Wide Association Study—GWAS

The genome wide association study was performed for milk production traits: (milk yield, milk protein concentration, and milk fat concentration) using eight different models (GLM, MLM, MLMM, FarmCPU, Super, ECMLM, CMLM, and Blink). The GWAS revealed an association of SNP rs110619097 located on BTA28 using GLM and Blink models with protein concentration (adjusted to 305 days) in the first lactation (Figure 1). We examined the Q–Q plots generated by the different models; Blink generated the line close to a 1:1 ratio with an upward deviated tail (Supplementary Figure S1). The SNP rs110619097 (28:g. 22479232T>C, UMD3.1.1) significantly (p < 0.01) associated with milk protein yield (adjusted to 305 days) is located on the BTA28 chromosome (Figure 2). The MAF for this SNP is 0.18 (p < 0.01). This marker is located within the intron of the CTNNA3 gene. The values of the milk protein concentration (adjusted to 305 days) in the first lactation for animals with different genotypes of the SNP rs110619097 were: CC: 3.263 ± 0.060 (10 animals), CT: 3.204 ± 0.109 (88 animals) and TT: 3.144 ± 0.134 (208 animals), undefined genotype: 3.320 (1 animal) (Figure 3). No other significant association of markers with investigated milk production traits (milk yield and milk fat concentration) was observed.

3.2. Runs of Homozygosity—ROH

After refining the SNP chip data, we used 52,934 markers to identify ROH segments in the genomes of the studied HF cows. We identified a total of 18,549 ROHs with 60.42 ROHs per animal ranging from 8 to 95. We classified ROHs into five different categories based on their length: 1–2 Mbp, 2–4 Mbp, 4–8 Mbp, 8–16 and >16 Mbps. Descriptive statistics for each length category are shown in Table 1. The majority of the detected ROHs were in the 1–2 Mbp length category (Table 1). The longest ROHs (>16 Mbp) were the least, with no ROHs in this category observed in 45% of animals (n = 139). The average length of the identified ROHs is 3.28 Mbp. The longest identified ROH is 65.31 Mbp long and is located on chromosome BTA7. The total length of all ROH segments per animal is 197.92 Mbp, where in the animal with the largest total ROH length (433.37 Mbp), 17% of the autosome is covered by an ROH. The value of genomic inbreeding (FROH) in that animal is 0.173, and in the least inbred animal, it is 0.005. Genomic inbreeding in the studied population is on average 0.079. For each length category, the inbreeding coefficient was calculated as the proportion of the autosomes of each individual that were covered by an ROH from that class (e.g., FROH > 16) (Figure 4).

3.3. Selection Signatures

In order to examine the effect of selection in the population of HF cows, we investigated the occurrence of ROH segments in the genomes of the examined individuals. The ROH incidence was not uniform across chromosomes (Figure 5). We identified ROH islands on chromosomes 1, 7, 10, 16, 22, and 26 (Supplementary Table S1). The highest ROH rate was observed on chromosome 1 (45% of animals) on markers AX-106742409 (rs41578805), AX-106740136 (rs109562914), AX-106719581 (rs41603780), AX-185115133 (rs523828967), and AX-124379394 (rs110792335).

4. Discussion

Various methods are being used to make breeding and selection programs more efficient [44]. Milk is a very important source of proteins and microelements and remains one of the most important components of the infant and adult diet. Dairy is a major dietary component of the human diet [45,46]. Having this in mind, the need for greater genetics of high-lactating cows is bigger than ever. Our study contributes to a better understanding of the genetics which underlie the milk productivity traits of cattle.
Even though GWASs have been conducted in cattle extensively in recent years, the available data are still insufficient to provide a complete explanation of the genetic mechanisms for milk production traits. As a result of numerous studies, hundreds of thousands of QTLs, which are responsible for the 680 traits [47], have been discovered on 30 chromosomes of cattle, and the list is still expanding. In our study, 334 HF cows with corresponding phenotypic data on milk yield, milk protein, and milk fat concentrations were genotyped and included in the GWAS. In this first of a kind study in Serbia, we found a significant correlation between the CTNNA3 gene on the BTA28 chromosome and milk protein concentrations in the first lactation (adjusted to 305 days). In the study conducted on American HF cows, Cole et al. [18], similarly to our investigation, discovered a significant association between the CTNNA3 gene and milk protein content (%), milk protein yield, milk fat content (%), milk fat yield, and milk yield. Using the GenABEL library for R software, v. 4.3.2, Dadousis et al. [48] observed the association of the same gene with other parameters that are important for dairy cattle (%CYSOLIDS, RECSOLIDS, and RECenergy—parameters related to curd dry matter). In general, many statistical methods and software packages are used for GWASs to improve computer efficiency, the power of statistical processing, and control of false positive results. In our research, we used GAPIT Version 3 software, which provided us with the possibility of comparing the results of association analyses using various statistical models (GLM, MLM, MLMM, FarmCPU, Super, ECMLM, CMLM, and Blink). In our study, the CTNNA3 gene exerted an association with milk protein concentration when each of these statistical models was applied, indicating the potential of this marker in future dairy cattle breeding and selection programs. It is interesting to mention that we did not observe any significant signals from some of the most commonly discussed genes (e.g., DGAT, CSN). This could be due to the specificity of our study, such as a relatively small number of samples obtained from only one farm.
The level of inbreeding in a population is an important parameter for monitoring its genetic diversity and its management. High levels of inbreeding cause inbreeding depression and should be avoided in farm animals [49]. Moreover, increasing inbreeding is associated with negative effects on the production and reproductive characteristics of dairy cows [50,51]. Inbreeding is traditionally measured using pedigree data. However, the pedigree-based inbreeding coefficient has certain limitations. Inbreeding is a widely accepted concept for characterizing the evolution, diversity, and general structure of a population, and it can be studied at different levels, starting from the individual, herd, and up to the population level [48,52] or even within a consortium, such as an artificial insemination center or dairy farm [53]. Furthermore, controlling the degree of inbreeding in cattle populations is particularly important, as only a small number of animals from the entire population are used for breeding. This knowledge allows for better management of animal diversity and genetics, controlling inbreeding depression and, in smaller populations, it can be used for conservation purposes. To calculate genomic inbreeding in our study, we used genomic information obtained by means of a medium-density SNP chip (63,648 markers). With the SNP data that are obtained, the expected inbreeding values are replaced by the resulting homozygosity measurements, which are thought to be a more accurate way to assess inbreeding and better reflect the level of homozygosity. In our study, the least detected ROHs were the longest (>16 Mbp—indicating recent inbreeding), and 45% of the animals had no ROH > 16 Mb at all, which is similar to the results obtained by Doekes et al. [54], who in their research did not observe ROH > 16 Mb in 26% of the examined cows. Our findings are in line with previous studies [41,55,56] in which short ROHs (<2 Mbp) were observed more frequently than long ROHs. Our results of mean ROH counts per animal, observed across all categories (1–2 Mbp: 30.140; 2–4 Mbp: 17.505; 4–8 Mbp: 8.679; 8–16 Mbp: 3.992; >16 Mbp: 1.726, respectively), are in accordance with the results of Marras et al. [55], who determined the following values: 1–2 Mbp: 46.5; 2–4 Mbps: 17.0; 4–8 Mbps: 9.7; 8–16 Mbps: 5.9; >16 Mbps: 3.0. In our study, the average genomic inbreeding value (FROH) was calculated using SVS Golden Helix software v8.9.0 and was 0.079, which is slightly lower than the values obtained by Marras et al. [55] in a population of Italian HF cattle (0.116), while Makanjuola et al. [57] obtained FROH values of 0.136 and 0.156, respectively, in a population of Canadian HF cattle, using SNP1101 version 1.0 and PLINK v1.90 software. Our results suggest the presence of a low level of recent and total inbreeding in the investigated population of Serbian HF cattle.
Identifying recent signatures of positive selection in domestic animals can provide information on the genomic regions affected by artificial and natural selection and thereby help identify beneficial mutations and underlying biological pathways for economically important traits [39]. There are several different approaches to detecting selection signatures, and we decided to utilize ROH analyses. ROHs are a state-of-the-art method for analyzing inbreeding in animal populations. Furthermore, ROHs are suitable for revealing selection signatures via ROH island identification. ROH islands can be defined as genomic regions with reduced genetic diversity and, consequently, high homozygosity around a selected locus that could contain targeted regions of positive selection, which are under strong selective pressure [43]. In our research, ROH islands were observed on chromosomes 1, 7, 10, 16, 22, and 26, while chromosome 1 had the highest rate of detected ROHs. ROH islands on chromosome BTA1 were also identified in other cattle populations by Purfield et al. [31], Mastrangelo et al. [58], and Gurgul et al. [59]. Our results indicate the existence of a 2 Mbp (starting at 83 Mbp and ending at 85 Mbp) ROH island on BTA1. In 45% of the animals from the studied population, selection signatures appear in the following genes: PARL, YEATS2, and KLHL24, while Serão et al. [60] observed an association of this region with residual feed intake (RFI) in the Angus and Simmental cattle breeds, which leads us to conclusion that the discovered ROH islands may be related to other productive/reproductive traits of cattle.
To summarize, by utilizing a GWAS, a significant association (p < 0.01) was observed between the rs110619097 polymorphism located in the intron of the CTNNA3 gene of the BTA28 chromosome and milk protein concentration in the first lactation (adjusted to 305 days). The results of the ROH segment analyses revealed a low level of recent and total genomic inbreeding in the examined population of HF cattle. The detected ROH islands imply their link with the productive traits of dairy cows and present the regions of the genome that are under selection pressure for other economically important traits.

5. Conclusions

This study presents a pilot study in the field of the SNP genotyping of animals in Serbia. All of the obtained data regarding the GWAS and ROH analyses provide useful information for a better understanding of the genetic architecture of milk production traits and genomic inbreeding in dairy cows and present a good base for many more studies on Serbian HF cattle. Furthermore, we are keen to continue our research on larger number of Serbian HF cattle, expanding it to larger number of farms and looking for other important traits of dairy cows, such as health and reproductive traits.

Supplementary Materials

Please add: The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/ani14050669/s1, Figure S1: Q-Q plot of GWAS for milk protein concentration in the first lactation (reduced to 305 days). Table S1: Detected selection signatures and the association found by other authors.

Author Contributions

Conceptualization, M.R., M.Z., U.G., M.M. and Z.S.; data curation, M.R., M.Z., U.G., J.B., M.M. and Z.S.; formal analysis, M.R., U.G. and J.B.; methodology, M.R., M.Z., U.G., J.S., M.M. and Z.S.; resources, Z.S.; software, M.Z.; supervision, Z.S.; writing—original draft, M.R. and M.Z.; writing—review and editing, U.G., J.S., M.M. and Z.S. All authors have read and agreed to the published version of the manuscript.

Funding

This study was supported and financed by: (i) the International Atomic Energy Agency (IAEA) through the Research Project entitled “Animal Identification, Pedigree, Exterior and Performance Data Recording in Selected Holstein-Friesian (HF) Cattle Population in Serbia Used for Future Genetic Selection under AI Programmes” (Research Contract 20774), which forms part of the IAEA Coordinated Research Project D31028, entitled “Application of Nuclear and Genomic Tools to Enable for the Selection of Animals with Enhanced Productivity Traits” and (ii) the Ministry of Science, Technological Development and Innovation of the Republic of Serbia (Contract No. 451-03-47/2023-01/200143). Both projects are led by Zoran Stanimirovic.

Institutional Review Board Statement

This research was conducted in accordance with all relevant national regulations and institutional policies for the care and use of animals—Ministry of Agriculture, Forestry and Water Management of the Republic of Serbia, Decision No: 323-07-06299/2018-05, issued on 23 July 2018.

Informed Consent Statement

The informed consent was obtained from the AI Dahra Serbia farm units.

Data Availability Statement

The raw data supporting the conclusions of this article will be made available by the authors on request.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Ristanic, M.; Stanisic, L.; Maletic, M.; Glavinic, U.; Draskovic, V.; Aleksic, N.; Stanimirovic, Z. Bovine foetal sex determination—Different DNA extraction and amplification approaches for efficient livestock production. Reprod. Domest. Anim. 2018, 53, 947–954. [Google Scholar] [CrossRef]
  2. Ristanić, M.; Glavinić, U.; Vejnović, B.; Maletić, M.; Kirovski, D.; Teodorović, V.; Stanimirović, Z. Beta-casein gene polymorphism in Serbian Holstein-Friesian cows and its relationship with milk production traits. Acta Vet. 2020, 70, 497–510. [Google Scholar]
  3. Pedrosa, V.B.; Schenkel, F.S.; Chen, S.Y.; Oliveira, H.R.; Casey, T.M.; Melka, M.G.; Brito, L.F. Genomewide association analyses of lactation persistency and milk production traits in Holstein cattle based on imputed whole-genome sequence data. Genes 2021, 12, 1830. [Google Scholar] [CrossRef] [PubMed]
  4. Bekele, R.; Taye, M.; Abebe, G.; Meseret, S. Genomic regions and candidate genes associated with milk production traits in Holstein and its crossbred cattle: A review. Int. J. Genom. 2023, 2023, 8497453. [Google Scholar] [CrossRef]
  5. Zhang, H.; Wang, Z.; Wang, S.; Li, H. Progress of genome wide association study in domestic animals. J. Anim. Sci. Biotechnol. 2012, 3, 26. [Google Scholar] [CrossRef] [PubMed]
  6. Saowaphak, P.; Duangjinda, M.; Plaengkaeo, S.; Suwannasing, R.; Boonkum, W. Genetic correlation and genome-wide association study (GWAS) of the length of productive life, days open, and 305-days milk yield in crossbred Holstein dairy cattle. Genet. Mol. Res. 2017, 16, 1–11. [Google Scholar] [CrossRef]
  7. Sorbolini, S.; Bongiorni, S.; Cellesi, M.; Gaspa, G.; Dimauro, C.; Valentini, A.; Macciotta, N.P.P. Genome wide association study on beef production traits in Marchigiana cattle breed. J. Anim. Breed. Genet. 2017, 134, 43–48. [Google Scholar] [CrossRef]
  8. Sermyagin, A.A.; Gladyr, E.A.; Plemyashov, K.V.; Kudinov, A.A.; Dotsev, A.V.; Deniskova, T.E.; Zinovieva, N.A. Genome-wide association studies for milk production traits in Russian population of Holstein and black-and-white cattle. In Proceedings of the Scientific-Practical Conference “Research and Development—2016”, Moscow, Russia, 14–15 December 2016; Springer: Berlin/Heidelberg, Germany, 2018; pp. 591–599. [Google Scholar]
  9. Zhou, J.; Liu, L.; Chen, C.J.; Zhang, M.; Lu, X.; Zhang, Z.; Huang, X.; Shi, Y. Genome-wide association study of milk and reproductive traits in dual-purpose Xinjiang Brown cattle. BMC Genom. 2019, 20, 827. [Google Scholar] [CrossRef]
  10. Paiva, J.T.; Peixoto, M.G.C.D.; Bruneli, F.A.T.; Alvarenga, A.B.; Oliveira, H.R.; Silva, A.A.; Silva, D.A.; Veroneze, R.; Silva, F.F.; Lopes, P.S. Genetic parameters, genome-wide association and gene networks for milk and reproductive traits in Guzerá cattle. Livest. Sci. 2020, 242, 104273. [Google Scholar] [CrossRef]
  11. Pegolo, S.; Momen, M.; Morota, G.; Rosa, G.J.; Gianola, D.; Bittante, G.; Cecchinato, A. Structural equation modeling for investigating multi-trait genetic architecture of udder health in dairy cattle. Sci. Rep. 2020, 10, 7751. [Google Scholar] [CrossRef]
  12. Da Cruz, A.S.; Silva, D.C.; Minasi, L.B.; de Farias Teixeira, L.K.; Rodrigues, F.M.; da Silva, C.C.; do Carmo, A.S.; da Silva, M.V.G.B.; Utsunomiya, Y.T.; Garcia, J.F.; et al. Single-Nucleotide Polymorphism variations associated with specific genes putatively identified enhanced genetic predisposition for 305-day milk yield in the Girolando crossbreed. Front. Genet. 2021, 11, 573344. [Google Scholar] [CrossRef] [PubMed]
  13. Kim, S.; Lim, B.; Cho, J.; Lee, S.; Dang, C.G.; Jeon, J.H.; Kim, J.M.; Lee, J. Genome-wide identification of candidate genes for milk production traits in Korean Holstein Cattle. Animals 2021, 11, 1392. [Google Scholar] [CrossRef] [PubMed]
  14. Lu, X.; Arbab, A.A.I.; Abdalla, I.M.; Liu, D.; Zhang, Z.; Xu, T.; Su, G.; Yang, Z. Genetic Parameter Estimation and Genome-Wide Association Study-Based Loci Identification of Milk-Related Traits in Chinese Holstein. Front. Genet. 2022, 12, 799664. [Google Scholar] [CrossRef]
  15. de Souza Fonseca, P.A.; Caldwell, T.; Mandell, I.; Wood, K.; Cánovas, A. Genome-wide association study for meat tenderness in beef cattle identifies patterns of the genetic contribution in different post-mortem stages. Meat Sci. 2022, 186, 108733. [Google Scholar] [CrossRef] [PubMed]
  16. Chen, S.Y.; Gloria, L.S.; Pedrosa, V.B.; Doucette, J.; Boerman, J.P.; Brito, L.F. Unravelling the genomic background of resilience based on variability in milk yield and milk production levels in North American Holstein cattle through GWAS and Mendelian randomization analyses. J. Dairy Sci. 2024, 107, 1035–1053. [Google Scholar] [CrossRef] [PubMed]
  17. Haque, M.A.; Alam, M.Z.; Iqbal, A.; Lee, Y.M.; Dang, C.G.; Kim, J.J. Genome-Wide Association Studies for Body Conformation Traits in Korean Holstein Population. Animals 2023, 13, 2964. [Google Scholar] [CrossRef] [PubMed]
  18. Cole, J.B.; Wiggans, G.R.; Ma, L.; Sonstegard, T.S.; Lawlor, T.J.; Crooker, B.A.; Van Tassell, C.P.; Yang, J.; Wang, S.; Matukumalli, L.K.; et al. Genome-wide association analysis of thirty one production, health, reproduction and body conformation traits in contemporary US Holstein cows. BMC Genom. 2011, 12, 408. [Google Scholar] [CrossRef] [PubMed]
  19. Kolbehdari, D.; Wang, Z.; Grant, J.R.; Murdoch, B.; Prasad, A.; Xiu, Z.; Marques, E.; Stothard, P.; Moore, S.S. A whole genome scan to map QTL for milk production traits and somatic cell score in Canadian Holstein bulls. J. Anim. Breed. Genet. 2009, 126, 216–227. [Google Scholar] [CrossRef]
  20. Jiang, J.; Ma, L.; Prakapenka, D.; VanRaden, P.M.; Cole, J.B.; Da, Y. A large-scale genome-wide association study in US Holstein cattle. Front. Genet. 2019, 10, 412. [Google Scholar] [CrossRef]
  21. Liu, L.; Zhou, J.; Chen, C.J.; Zhang, J.; Wen, W.; Tian, J.; Zhang, Z.; Gu, Y. GWAS-based identification of new loci for milk yield, fat, and protein in Holstein cattle. Animals 2020, 10, 2048. [Google Scholar] [CrossRef]
  22. Schopen, G.; Visker, M.; Koks, P.; Mullaart, E.; van Arendonk, J.; Bovenhuis, H. Whole-genome association study for milk protein composition in dairy cattle. J. Dairy Sci. 2011, 94, 3148–3158. [Google Scholar] [CrossRef] [PubMed]
  23. Fontanesi, L.; Calò, D.G.; Galimberti, G.; Negrini, R.; Marino, R.; Nardone, A.; Ajmone-Marsan, P.; Russo, V. A candidate gene association study for nine economically important traits in Italian Holstein cattle. Anim. Genet. 2014, 45, 576–580. [Google Scholar] [CrossRef]
  24. Jiang, J.; Liu, L.; Gao, Y.; Shi, L.; Li, Y.; Liang, W.; Sun, D. Determination of genetic associations between indels in 11 candidate genes and milk composition traits in Chinese Holstein population. BMC Genet. 2019, 20, 48. [Google Scholar] [CrossRef] [PubMed]
  25. Murgiano, L.; Jagannathan, V.; Benazzi, C.; Bolcato, M.; Brunetti, B.; Muscatello, L.V.; Dittmer, K.; Piffer, C.; Gentile, A.; Drögemüller, C. Deletion in the EVC2 gene causes chondrodysplastic dwarfism in Tyrolean Grey cattle. PLoS ONE 2014, 9, e94861. [Google Scholar] [CrossRef] [PubMed]
  26. Tiezzi, F.; Parker-Gaddis, K.L.; Cole, J.B.; Clay, J.S.; Maltecca, C. A genome-wide association study for clinical mastitis in first parity US Holstein cows using single-step approach and genomic matrix re-weighting procedure. PLoS ONE 2015, 10, e0114919. [Google Scholar] [CrossRef]
  27. Brown, A.; Ojango, J.; Gibson, J.; Coffey, M.; Okeyo, M.; Mrode, R. Genomic selection in a crossbred cattle population using data from the dairy genetics East Africa project. J. Dairy Sci. 2016, 99, 7308–7312. [Google Scholar] [CrossRef]
  28. Boison, S.A.; Utsunomiya, A.T.H.; Santos, D.J.A.; Neves, H.H.R.; Carvalheiro, R.; Mészáros, G.; Utsunomiya, Y.T.; do Carmo, A.S.; Verneque, R.S.; Machado, M.A.; et al. Accuracy of genomic predictions in Gyr (Bos indicus) dairy cattle. J. Dairy Sci. 2017, 100, 5479–5490. [Google Scholar] [CrossRef]
  29. Meuwissen, T.; Goddard, M. The use of family relationships and linkage disequilibrium to impute phase and missing genotypes in up to whole-genome sequence density genotypic data. Genetics 2010, 185, 1441–1449. [Google Scholar] [CrossRef]
  30. VanRaden, P.M. Efficient methods to compute genomic predictions. J. Dairy Sci. 2008, 91, 4414–4423. [Google Scholar] [CrossRef]
  31. Purfield, D.C.; Berry, D.P.; McParland, S.; Bradley, D.G. Runs of homozygosity and population history in cattle. BMC Genet. 2012, 13, 70. [Google Scholar] [CrossRef]
  32. Ferenčaković, M.; Hamzić, E.; Gredler, B.; Solberg, T.R.; Klemetsdal, G.; Curik, I.; Sölkner, J. Estimates of autozygosity derived from runs of homozygosity: Empirical evidence from selected cattle populations. J. Anim. Breed. Genet. 2013, 130, 286–293. [Google Scholar] [CrossRef] [PubMed]
  33. Ferenčaković, M.; Sölkner, J.; Curik, I. Estimating autozygosity from high-throughput information: Effects of SNP density and genotyping errors. Genet. Sel. Evol. 2013, 45, 42. [Google Scholar] [CrossRef]
  34. Zhang, Q.; Calus, M.P.; Guldbrandtsen, B.; Lund, M.S.; Sahana, G. Estimation of inbreeding using pedigree, 50k SNP chip genotypes and full sequence data in three cattle breeds. BMC Genet. 2015, 16, 88. [Google Scholar] [CrossRef] [PubMed]
  35. Zhang, Q.; Guldbrandtsen, B.; Bosse, M.; Lund, M.S.; Sahana, G. Runs of homozygosity and distribution of functional variants in the cattle genome. BMC Genom. 2015, 16, 542. [Google Scholar] [CrossRef] [PubMed]
  36. Makanjuola, B.O.; Maltecca, C.; Miglior, F.; Marras, G.; Abdalla, E.A.; Schenkel, F.S.; Baes, C.F. Identification of unique ROH regions with unfavorable effects on production and fertility traits in Canadian Holsteins. Genet. Sel. Evol. 2021, 53, 68. [Google Scholar] [CrossRef] [PubMed]
  37. Rothammer, S.; Seichter, D.; Förster, M.; Medugorac, I. A genome-wide scan for signatures of differential artificial selection in ten cattle breeds. BMC Genom. 2013, 14, 908. [Google Scholar] [CrossRef] [PubMed]
  38. Gouveia, J.J.D.S.; Silva, M.V.G.B.D.; Paiva, S.R.; Oliveira, S.M.P.D. Identification of selection signatures in livestock species. Genet. Mol. Biol. 2014, 37, 330–342. [Google Scholar] [CrossRef]
  39. Ristanic, M. Comparative Analyses of β-Casein Gene Polymorphism (A1/A2 Genotype) and Its Effect on Qualitative Composition of Milk in Lactating and Autochthonous Cattle Breeds. Ph.D. Thesis, Faculty of Veterinary Medicine, University of Belgrade, Belgrade, Serbia, 2022. [Google Scholar]
  40. Wang, J.; Zhang, Z. GAPIT version 3: Boosting power and accuracy for genomic association and prediction. Genom. Proteom. Bioinform. 2021, 19, 629–640. [Google Scholar] [CrossRef]
  41. McQuillan, R.; Leutenegger, A.L.; Abdel-Rahman, R.; Franklin, C.S.; Pericic, M.; Barac-Lauc, L.; Smolej-Narancic, N.; Janicijevic, B.; Polasek, O.; Tenesa, A.; et al. Runs of homozygosity in European populations. Am. J. Hum. Genet. 2008, 83, 359–372. [Google Scholar] [CrossRef]
  42. R Core Team. R: A Language and Environment for Statistical Computing; The R Foundation: Vienna, Austria, 2013. [Google Scholar]
  43. Gorssen, W.; Meyermans, R.; Janssens, S.; Buys, N. A publicly available repository of ROH islands reveals signatures of selection in different livestock and pet species. Genet. Sel. Evol. 2021, 53, 2. [Google Scholar] [CrossRef]
  44. Sahana, G.; Cai, Z.; Sanchez, M.P.; Bouwman, A.C.; Boichard, D. Invited review: Good practices in genome-wide association studies to identify candidate sequence variants in dairy cattle. J. Dairy Sci. 2023, 106, 5218–5241. [Google Scholar] [CrossRef]
  45. Prylipko, T.; Kuzminska, I.; Sheludchenko, L.; Semenyshena, R.; Koval, T. Control of the Quality and Safety of Dairy Products in Ukraine: International and Legal Aspects. Eur. Food Feed Law Rev. 2023, 18, 22. [Google Scholar]
  46. Rahimzadeh Barzoki, H.; Faraji, H.; Beirami, S.; Keramati, F.Z.; Nayik, G.A.; Izadi Yazdanaabadi, Z.; Mozaffari Nejad, A.S. Seasonal Study of Aflatoxin M1 Contamination in Cow Milk on the Retail Dairy Market in Gorgan, Iran. Dairy 2023, 4, 571–580. [Google Scholar] [CrossRef]
  47. Cattle QTLdb. Available online: https://www.animalgenome.org/cgi-bin/QTLdb/BT/qdetails?QTL_ID=65988 (accessed on 10 December 2023).
  48. Dadousis, C.; Pegolo, S.; Rosa, G.J.M.; Gianola, D.; Bittante, G.; Cecchinato, A. Pathway-based genome-wide association analysis of milk coagulation properties, curd firmness, cheese yield, and curd nutrient recovery in dairy cattle. J. Dairy Sci. 2017, 100, 1223–1231. [Google Scholar] [CrossRef] [PubMed]
  49. Nishio, M.; Inoue, K.; Ogawa, S.; Ichinoseki, K.; Arakawa, A.; Fukuzawa, Y.; Okamura, T.; Kobayashi, E.; Taniguchi, M.; Oe, M.; et al. Comparing pedigree and genomic inbreeding coefficients, and inbreeding depression of reproductive traits in Japanese Black cattle. BMC Genom. 2023, 24, 376. [Google Scholar] [CrossRef]
  50. Bjelland, D.W.; Weigel, K.A.; Vukasinovic, N.; Nkrumah, J.D. Evaluation of inbreeding depression in Holstein cattle using whole-genome SNP markers and alternative measures of genomic inbreeding. J. Dairy Sci. 2013, 96, 4697–4706. [Google Scholar] [CrossRef]
  51. Lozada-Soto, E.A.; Maltecca, C.; Lu, D.; Miller, S.; Cole, J.B.; Tiezzi, F. Trends in genetic diversity and the effect of inbreeding in American Angus cattle under genomic selection. Genet. Sel. Evol. 2021, 53, 50. [Google Scholar] [CrossRef] [PubMed]
  52. Howard, J.T.; Pryce, J.E.; Baes, C.; Maltecca, C. Invited review: Inbreeding in the genomics era: Inbreeding, inbreeding depression, and management of genomic variability. J. Dairy Sci. 2017, 100, 6009–6024. [Google Scholar] [CrossRef] [PubMed]
  53. Ablondi, M.; Malacarne, M.; Cipolat-Gotet, C.; van Kaam, J.T.; Sabbioni, A.; Summer, A. Genome-wide scan reveals genetic divergence in Italian Holstein cows bred within PDO cheese production chains. Sci. Rep. 2021, 15, 12601. [Google Scholar] [CrossRef]
  54. Doekes, H.P.; Veerkamp, R.F.; Bijma, P.; de Jong, G.; Hiemstra, S.J.; Windig, J.J. Inbreeding depression due to recent and ancient inbreeding in Dutch Holstein–Friesian dairy cattle. Genet. Sel. Evol. 2019, 51, 54. [Google Scholar] [CrossRef]
  55. Marras, G.; Gaspa, G.; Sorbolini, S.; Dimauro, C.; Ajmone-Marsan, P.; Valentini, A.; Williams, J.L.; Macciotta, N.P. Analysis of runs of homozygosity and their relationship with inbreeding in five cattle breeds farmed in Italy. Anim. Genet. 2015, 46, 110–121. [Google Scholar] [CrossRef] [PubMed]
  56. Forutan, M.; Ansari Mahyari, S.; Baes, C.; Melzer, N.; Schenkel, F.S.; Sargolzaei, M. Inbreeding and runs of homozygosity before and after genomic selection in North American Holstein cattle. BMC Genom. 2018, 19, 98. [Google Scholar] [CrossRef] [PubMed]
  57. Makanjuola, B.O.; Miglior, F.; Abdalla, E.A.; Maltecca, C.; Schenkel, F.S.; Baes, C.F. Effect of genomic selection on rate of inbreeding and coancestry and effective population size of Holstein and Jersey cattle populations. J. Dairy Sci. 2020, 103, 5183–5199. [Google Scholar] [CrossRef] [PubMed]
  58. Mastrangelo, S.; Tolone, M.; Di Gerlando, R.; Fontanesi, L.; Sardina, M.T.; Portolano, B. Genomic inbreeding estimation in small populations: Evaluation of runs of homozygosity in three local dairy cattle breeds. Animal 2016, 10, 746–754. [Google Scholar] [CrossRef]
  59. Gurgul, A.; Szmatoła, T.; Topolski, P.; Jasielczuk, I.; Żukowski, K.; Bugno-Poniewierska, M. The use of runs of homozygosity for estimation of recent inbreeding in Holstein cattle. J. Appl. Genet. 2016, 57, 527–530. [Google Scholar] [CrossRef]
  60. Serão, N.V.; González-Peña, D.; Beever, J.E.; Faulkner, D.B.; Southey, B.R.; Rodriguez-Zas, S.L. Single nucleotide polymorphisms and haplotypes associated with feed efficiency in beef cattle. BMC Genet. 2013, 14, 94. [Google Scholar] [CrossRef]
Figure 1. The circular Manhattan plot shows the results of eight different models: (I) GLM, (II) MLM, (III) MLMM, (IV) FarmCPU, (V) Super, (VI) ECMLM, (VII) CMLM, and (VIII) Blink; GWAS analysis of milk protein concentration (adjusted to 305 days) in the first lactation. The overlap in results between two different models is shown as a dashed vertical line extending through the Manhattan plots for all models used.
Figure 1. The circular Manhattan plot shows the results of eight different models: (I) GLM, (II) MLM, (III) MLMM, (IV) FarmCPU, (V) Super, (VI) ECMLM, (VII) CMLM, and (VIII) Blink; GWAS analysis of milk protein concentration (adjusted to 305 days) in the first lactation. The overlap in results between two different models is shown as a dashed vertical line extending through the Manhattan plots for all models used.
Animals 14 00669 g001
Figure 2. Manhattan plot of the GWAS analysis for the milk protein concentration (adjusted to 305 days) in the first lactation. The plot represents the −log10 transformed p values for all analyzed SNPs, green line presents GWAS significance level.
Figure 2. Manhattan plot of the GWAS analysis for the milk protein concentration (adjusted to 305 days) in the first lactation. The plot represents the −log10 transformed p values for all analyzed SNPs, green line presents GWAS significance level.
Animals 14 00669 g002
Figure 3. Milk protein concentration (adjusted to 305 days) in the first lactation with three different genotypes in the rs110619097 (28:g. 22479232T> C). X axis: different genotypes; Y axis: milk protein concentration.
Figure 3. Milk protein concentration (adjusted to 305 days) in the first lactation with three different genotypes in the rs110619097 (28:g. 22479232T> C). X axis: different genotypes; Y axis: milk protein concentration.
Animals 14 00669 g003
Figure 4. Boxplot distribution of genomic inbreeding (FROH) obtained from different ROH length categories.
Figure 4. Boxplot distribution of genomic inbreeding (FROH) obtained from different ROH length categories.
Animals 14 00669 g004
Figure 5. Manhattan plot distribution of ROH islands on chromosomes 1, 7, 10, 16, 22, and 26. The X-axis represents the chromosomal distribution of ROHs across the genome; The Y-axis shows the frequency of overlapping ROHs that are common in the tested samples.
Figure 5. Manhattan plot distribution of ROH islands on chromosomes 1, 7, 10, 16, 22, and 26. The X-axis represents the chromosomal distribution of ROHs across the genome; The Y-axis shows the frequency of overlapping ROHs that are common in the tested samples.
Animals 14 00669 g005
Table 1. ROH length, ROH number, and genomic inbreeding (FROH) in different ROH length categories per animal.
Table 1. ROH length, ROH number, and genomic inbreeding (FROH) in different ROH length categories per animal.
ROH Length Category
Parameter1–2 Mbp2–4 Mbp4–8 Mbp8–16 Mbp>16 Mbp
Total length of ROH per animalAverage Value43.84848.11948.14544.03839.332
Standard Deviation9.97915.50420.13825.65524.867
Median43.81748.55347.65139.14428.769
Minimum7.4162.9274.2448.01316.282
Maximum66.53893.096119.858144.577141.065
Number of ROH per animalAverage Value30.14017.5058.6793.9921.726
Standard Deviation6.7235.4883.5502.2860.933
Median30.00018.0008.0004.0001.000
Minimum6.0001.0001.0001.0001.000
Maximum45.00035.00021.00012.0005.000
Genomic inbreeding FROHAverage Value0.0170.0190.0190.0180.016
Standard Deviation0.0040.0060.0080.0100.010
Median0.0170.0190.0190.0160.011
Minimum0.0030.0010.0020.0030.006
Maximum0.0260.0370.0480.0580.056
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Ristanic, M.; Zorc, M.; Glavinic, U.; Stevanovic, J.; Blagojevic, J.; Maletic, M.; Stanimirovic, Z. Genome-Wide Analysis of Milk Production Traits and Selection Signatures in Serbian Holstein-Friesian Cattle. Animals 2024, 14, 669. https://doi.org/10.3390/ani14050669

AMA Style

Ristanic M, Zorc M, Glavinic U, Stevanovic J, Blagojevic J, Maletic M, Stanimirovic Z. Genome-Wide Analysis of Milk Production Traits and Selection Signatures in Serbian Holstein-Friesian Cattle. Animals. 2024; 14(5):669. https://doi.org/10.3390/ani14050669

Chicago/Turabian Style

Ristanic, Marko, Minja Zorc, Uros Glavinic, Jevrosima Stevanovic, Jovan Blagojevic, Milan Maletic, and Zoran Stanimirovic. 2024. "Genome-Wide Analysis of Milk Production Traits and Selection Signatures in Serbian Holstein-Friesian Cattle" Animals 14, no. 5: 669. https://doi.org/10.3390/ani14050669

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop