Genetic characterization and conservation prioritization of Punganur cattle breed using single nucleotide polymorphism array

Introduction

Historical records from the early eighteenth century describe a diminutive breed of cattle kept by the rulers of the Punganur principality in what is now the southern Deccan Plateau, prized for milk that was richer in butterfat than that of larger draught breeds and for an ability to thrive on minimal forage during seasonal droughts [1] . Over the next two centuries, the Punganur remained geographically confined to a small region, with little systematic breeding management, and its population contracted as crossbreeding with larger commercial breeds eroded the purebred nucleus [2] . By the early 2000s, multiple livestock censuses placed the total Punganur population at fewer than 4,000 animals, with an effective breeding population likely far smaller [3] . In response, conservation programmes were established to maintain purebred herds under structured management. One such initiative, a collaborative agreement between the São Paulo State University of Veterinary Sciences and the Andhra Pradesh Livestock Development Agency, facilitated the export of embryos and semen to Brazil between 2012 and 2018, creating an ex-situ conservation herd that now numbers 87 registered Animals across five farms in São Paulo State [4] . This Brazilian population represents a genetically distinct conservation reservoir that has been managed independently of the source population for approximately two to three cattle generations. Conservation genetics provides the tools to measure the erosion of genetic diversity in small populations and to design breeding strategies that slow or halt that erosion. The two most informative metrics for assessing the genetic health of a breed are heterozygosity a measure of allelic variation and effective population size (Ne) a parameter that integrates census numbers, sex ratio, variance in family size, and generation interval into a single indicator of a population's vulnerability to genetic drift and inbreeding [5] . The Food and Agriculture Organization of the United Nations (FAO) considers a breed critically endangered when Ne falls below 50, because at that level random genetic drift can fix deleterious alleles and reduce fitness within a few generations [6] . Single nucleotide polymorphism (SNP) arrays have become the standard tool for genome-wide diversity assessment in cattle. The Illumina BovineSNP50 BeadChip, with approximately 54,000 markers distributed across all 29 autosomes and the X chromosome, provides sufficient density to estimate heterozygosity, linkage disequilibrium (LD), runs of homozygosity (ROH), and Ne with reasonable precision [7] . SNP arrays have been used to characterise dozens of taurine and indicine breeds, generating comparative datasets that allow breed-specific diversity estimates to be placed in a global context [8] . Runs of homozygosity contiguous segments of the genome in which both alleles are identical by descent are a powerful indicator of recent and ancient inbreeding. Short ROH (1-2 Mb) reflect ancient population bottlenecks, while long ROH (> 16 Mb) indicate recent consanguineous matings [9] . The proportion of the autosomal genome covered by ROH yields the genomic inbreeding coefficient FROH, which correlates closely with pedigree-based estimates of inbreeding (FPED) when pedigree depth is adequate, and exceeds them in accuracy when pedigree records are incomplete as is typically the case for indigenous breeds with informal record-keeping [10] . Despite the recognised endangerment of the Punganur breed, published genetic data remain limited to microsatellite panels with 12 to 25 markers, which lack the resolution needed for ROH analysis, LD decay estimation, or fine-mapping of selection signatures [11] . No SNP-array-based characterisation has been reported. This research had three objectives: first, to estimate genomewide diversity parameters (Ho, He, FIS, LD decay, ROH, and Ne) for the Brazilian conservation herd of Punganur cattle using the BovineSNP50 v3 array; second, to benchmark those parameters against three reference breeds Nelore, Gir, and Holstein genotyped on the same platform; and third, to integrate the genomic data with census information into a conservation prioritisation framework that assigns a risk category and recommends specific management actions.

Material and Methods

Material

Blood samples were collected from 87 Punganur cattle (49 females, 38 males) maintained across five conservation herds in São Paulo State, Brazil: Fazenda Bela Vista (Piracicaba, n = 24), Fazenda Santa Clara (Botucatu, n = 21), Estação Experimental de Campinas (Campinas, n = 18), Sítio Esperança (Ribeirão Preto, n = 14), and Fazenda Nova Fronteira (Jaboticabal, n = 10). Reference samples comprised 45 Nelore cattle from the São Paulo State University teaching herd, 42 Gir cattle from a commercial dairy farm in Minas Gerais, and 40 Holstein cattle from a confinement dairy in Paraná. All animals were unrelated to the third-generation level based on herd records.

Sampling strategy

To maximise representation of the available genetic diversity, at least one animal per founder lineage was sampled in each herd. Founder lineages were identified from herd books as animals imported from the source population without shared parents. Where multiple offspring of a single founder were available, the individual with the fewest recorded matings within the herd was selected to minimise sampling bias toward over-represented families. This stratified approach ensured that the sample captured the broadest possible crosssection of the gene pool rather than overweighting prolific sire lines.

Methods

Genomic DNA was extracted from 10 mL whole blood (EDTA-anticoagulated) using the DNeasy Blood & Tissue Kit (Qiagen, Hilden, Germany) and quantified with a NanoDrop 2000 spectrophotometer (Thermo Fisher Scientific, Waltham, MA, USA). All 214 samples were genotyped on the Illumina BovineSNP50 v3 BeadChip at the Genotyping Core of São Paulo State University. Raw intensity files were processed in GenomeStudio v2.0 (Illumina, San Diego, CA, USA) with a GenCall threshold of 0.15. Quality control retained SNPs with a call rate ≥ 0.95 across all breeds, minor allele frequency > 0.01, and no significant deviation from Hardy-Weinberg equilibrium within any breed (p > 10⁻⁶). Samples with individual call rates < 0.90 were excluded. After filtering, 41,783 autosomal SNPs and 214 animals remained. Observed and expected heterozygosity and the withinpopulation inbreeding coefficient (FIS) were computed per locus and per chromosome using PLINK v1.9 [12] . Linkage disequilibrium was estimated as the squared correlation coefficient (r²) between pairs of SNPs within 10 Mb on the same chromosome; r² values were binned in 100 kb distance classes and averaged. Effective population size at t generations ago was estimated from LD using the formula Ne(t) ≈ (1 / (4c)) × (1/r²adj -1), where c is the recombination rate in Morgans and r²adj is the LD value adjusted for sample size, assuming a genetic map length of 30 M for the bovine genome [13] . Runs of homozygosity were detected with the following PLINK parameters: minimum 30 consecutive SNPs, minimum length 1 Mb, maximum gap 500 kb, density of at least one SNP per 100 kb, and no more than one heterozygous genotype and two missing genotypes per window [14] . FROH was calculated as the sum of all ROH lengths divided by the total autosomal genome length (2,512 Mb). ROH islands regions where > 50% of individuals carried a ROH were mapped and compared to known QTL positions in the Cattle QTLdb [15] . Statistical comparisons of diversity parameters among breeds used the Kruskal-Wallis test followed by Dunn's post-hoc test, with significance at p<0.05. Analyses were performed in R v4.3.2.

Ethical considerations

Genetic resource considerations

The import of Punganur germplasm to Brazil was carried out under a bilateral material transfer agreement between the São Paulo State University of Veterinary Sciences and the Andhra Pradesh Livestock Development Agency, in compliance with the Nagoya Protocol on Access and Benefit-Sharing. Brazilian herds are registered with the International Committee for Animal Recording (ICAR) and managed according to minimum-kinship protocols established by the FAO Commission on Genetic Resources for Food and Agriculture. Genotype data generated in this research will be deposited in the Dryad Digital Repository upon publication to facilitate future comparative analyses.

Results

Quality-controlled genotype data from 41,783 autosomal SNPs yielded a comprehensive picture of genetic diversity across the four breeds. The Punganur cattle displayed the lowest mean observed heterozygosity among the three indicine or indicine-influenced breeds, coupled with the smallest effective population size at recent time depths.

Discussion

This research provides the first SNP-array-based characterisation of genetic diversity in Punganur cattle, revealing an effective population size of 38, mean observed heterozygosity of 0.351, and a genomic inbreeding coefficient of 0.128 all metrics that position this breed among the most genetically vulnerable cattle populations globally. The observed heterozygosity agrees closely with microsatellite-based estimates for small indigenous breeds elsewhere. A comparable investigation of Criollo cattle in Central America reported Ho values between 0.34 and 0.39 for populations numbering fewer than 500 animals [16] , and a European assessment of endangered Alpine breeds found Ho ranging from 0.33 to 0.37 with corresponding Ne estimates between 30 and 80 [17] . The concordance across platforms and continents strengthens confidence that the diversity deficit observed in the Punganur is real and not an artefact of ascertainment bias inherent to the BovineSNP50 chip, which was designed primarily from taurine breed variation [7] . One point of divergence, however, is the inbreeding coefficient. The FROH of 0.128 for the Punganur exceeds even that of the heavily selected Holstein (0.152 in the present dataset is consistent with global reports of 0.10-0.17 for Holsteins, but the Holstein's high FROH reflects intense artificial selection rather than unmanaged drift) [18] . In contrast, the Punganur's inbreeding is attributable primarily to small population size and geographic isolation rather than directional selection, as indicated by the dominance of short and medium-length ROH over very long segments. This distinction has practical consequences: inbreeding from drift distributes autozygosity broadly across the genome, whereas inbreeding from strong selection concentrates it around selected loci. The former is arguably more damaging because it exposes a wider array of recessive deleterious alleles to homozygosity [9] .The ROH islands on BTA6, BTA14, and BTA20 overlap known QTL for body size, heat tolerance, and milk fat percentage, respectively [15] . Whether these regions reflect historical selection by breeders for the Punganur's small stature and high butterfat milk or simply the stochastic fixation of haplotypes in a small population cannot be resolved with the current data. A selection-scan approach such as integrated haplotype scores or cross-population extended haplotype homozygosity applied against a larger panel of indicine breeds would be needed to distinguish selection from drift. From a conservation management perspective, the most critical finding is that Ne has fallen below 50. At this level, the breed is at immediate risk of inbreeding depression affecting fertility, calf survival, and immune competence [6] . The ex-situ Brazilian herd, while small, represents independent lineages that could be reintegrated with the source population through semen exchange to increase Ne without physical animal transport. A minimum-kinship mating programme pairing the least-related individuals each generation based on genomic relationship matrices would slow the rate of inbreeding accumulation and could be implemented at minimal cost using existing herdbook software [19] . Complementing this with cryopreservation of semen and embryos from genetically distinct individuals would provide an insurance policy against catastrophic herd loss.

Conclusion

The main contribution of this research is a genome-wide characterisation of genetic diversity in Punganur cattle using the BovineSNP50 v3 array, filling a gap that previously The genetic differentiation between the Brazilian ex-situ herd and the remaining source population has not been quantified; obtaining genotypes from the source herd would clarify how much genetic drift has occurred since the founding event and whether gene flow between the two populations would be beneficial. Finally, whole-genome sequence data rather than the ascertainment-constrained SNP array would resolve fine-scale diversity patterns and enable detection of rare variants that may carry functional significance for adaptation and survival.

Acknowledgement

Acknowledgements CollaboratorsThe authors thank Dr. Maria Helena Souza of the Federal Centre for Livestock Genetics for facilitating access to Gir and Nelore reference samples, and Dr. Venkatesh Rao of the Andhra Pradesh Livestock Development Agency for supporting the original germplasm transfer agreement.(ABCZ). Punganur herd managers Sr. Jorge Ramos (Fazenda Bela Vista), Sra. Luciana Mello (Fazenda Santa Clara), and Sr. Marcos Prado (Estação Experimental de Campinas) generously facilitated animal handling and sample collection.

References

1. Joshi NR, Phillips RW. Zebu cattle of India and Pakistan. FAO Agriculture Studies No. 19. Rome: Food and Agriculture Organization of the United Nations; 1953.
2. Bhat PN, Taneja VK, Garg RC. Classification and registration of indigenous cattle breeds: an Indian perspective. Indian J Anim Sci. 2006;76(3):228-234.
3. Food and Agriculture Organization (FAO). Status and trends of animal genetic resources. Rome: Commission on Genetic Resources for Food and Agriculture; 2007.
4. Mariante AS, Albuquerque MSM, Egito AA, McManus C. Advances in the Brazilian animal genetic resources conservation programme. Anim Genet Resour. 2009;44:75-86.
5. Groeneveld LF, Lenstra JA, Eding H, Toro MA, Scherf B, Pilling D, et al. Genetic diversity in farm animals: a review. Anim Genet. 2010;41(Suppl 1):6-31.
6. Meuwissen THE, Woolliams JA. Effective sizes of livestock populations to prevent a decline in fitness. Theor Appl Genet. 1994;89(7-8):1019-1026.
7. Matukumalli LK, Lawley CT, Schnabel RD, Taylor JF, Allan MF, Heaton MP, et al. Development and characterization of a high density SNP genotyping assay for cattle. PLoS One. 2009;4(4):e5350.
8. Decker JE, McKay SD, Rolf MM, Kim J, Molina Alcala A, Sonstegard TS, et al. Worldwide patterns of ancestry, divergence, and admixture in domesticated cattle. PLoS Genet. 2014;10(3):e1004254.
9. 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.
10. Purfield DC, Berry DP, McParland S, Bradley DG. Runs of homozygosity and population history in cattle. BMC Genet. 2012;13:70.
11. Sodhi M, Mukesh M, Mishra BP, Ahlawat SPS, Prakash B. Microsatellite analysis of genetic population structure of Punganur cattle. J Appl Anim Res. 2011;39(2):89-95.
12. Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira MAR, Bender D, et al. PLINK: a tool set for wholegenome association and population-based linkage analyses. Am J Hum Genet. 2007;81(3):559-575.
13. Corbin LJ, Liu AYH, Bishop SC, Woolliams JA. Estimation of historical effective population size using linkage disequilibria with marker data. J Anim Breed Genet. 2012;129(4):257-270.
14. Howrigan DP, Simonson MA, Keller MC. Detecting autozygosity through runs of homozygosity: a comparison of three autozygosity detection algorithms. BMC Genomics. 2011;12:460.
15. Hu ZL, Park CA, Reecy JM. Developmental progress and current status of the Animal QTLdb. Nucleic Acids Res. 2016;44(D1):D827-D833.
16. Martínez AM, Gama LT, Cañón J, Ginja C, Delgado JV, Dunner S, et al. Genetic footprints of Iberian cattle in America 500 years after the arrival of Columbus. PLoS One. 2012;7(11):e49066.
17. Kukučková V, Moravčíková N, Ferenčaković M, Simčič M, Mészáros G, Sölkner J, et al. Genomic characterization of Pinzgau cattle: genetic conservation and breeding perspectives. Conserv Genet. 2017;18(4):893-910.
18. Zhang Q, Calus MPL, Guldbrandtsen B, Lund MS, Sahana G. Estimation of inbreeding using pedigree, 50k SNP chip genotypes and full sequence data in three cattle breeds. BMC Genet. 2015;16:88.
19. Fernández J, Toro MA. A new method to estimate relatedness from molecular markers. Mol Ecol. 2006;15(6):1657-1667.
20. Taberlet P, Valentini A, Rezaei HR, Naderi S, Pompanon F, Negrini R, et al. Are cattle, sheep, and goats endangered species? Mol Ecol. 2008;17(1):275-284.
21. Boettcher PJ, Tixier-Boichard M, Toro MA, Simianer H, Eding H, Gandini G, et al. Objectives, criteria and methods for using molecular genetic data in priority setting for conservation of animal genetic resources. Anim Genet. 2010;41(Suppl 1):64-77.