Front. Anim. Sci.Frontiers in Animal ScienceFront. Anim. Sci.2673-6225Frontiers Media S.A.10.3389/fanim.2023.1249470Animal ScienceOriginal ResearchEstimation of genetic parameters for parasite resistance and genome-wide identification of runs of homozygosity islands in Florida Cracker sheepHidalgoJorge1†Estrada-ReyesZaira M.23*†OgunadeIbukun M.4Pech-CervantesAndres A.5TerrillThomas6IdowuModoluwamu D.4TaiwoGodstime41Department of Animal and Dairy Science, University of Georgia, Athens, GA, United States2Department of Animal Science, North Carolina Agricultural and Technical State University, Greensboro, NC, United States3Department of Animal Sciences, University of Florida, Gainesville, FL, United States4Division of Animal and Nutritional Science, West Virginia University, Morgantown, WV, United States5Department of Agriculture, Food and Resource Sciences, University of Maryland Eastern Shore, Princess Anne, MD, United States6College of Agricultural, Family Sciences, and Technology, Fort Valley State University, Fort Valley, GA, United States
Edited by: Beatriz Gutiérrez-Gil, University of León, Spain
Reviewed by: Shamik Polley, West Bengal University of Animal and Fishery Sciences, India
Onur Yilmaz, Adnan Menderes University, Türkiye
Valentina Riggio, University of Edinburgh, United Kingdom
*Correspondence: Zaira M. Estrada-Reyes, zmestradareyes@ncat.edu
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
In this study, we estimated genetic parameters for parasite resistance traits and the distribution of runs of homozygosity islands in Florida Cracker sheep. The dataset contained 365 animals with phenotypic records at 38 days post-infection for fecal egg count (FEC), packed cell volume (PCV), FAMACHA score (FAM), and body condition score (BCS). The pedigree file contained 695 animals born between 2016 and 2020 and included 279 individuals with genotypes. Genetic parameters were estimated using a multi-trait model with a Bayesian implementation via Gibbs sampling in the GIBBS3F90 program. Heritability was 0.33 ± 0.09 for FEC, 0.31 ± 0.10 for FAM, 0.22 ± 0.09 for PCV, and 0.19 ± 0.07 for BCS. The genetic correlation between FEC and FAM was 0.51 ± 0.21; the remaining genetic correlations had large posterior standard deviations and yielded 95% posterior intervals including zero or with values out of the parameter space because of our small dataset. Analysis of the distribution of runs of homozygosity islands revealed 113 hot spots with annotated genes related to immune response and parasite resistance traits. Our results suggest that the genetic selection for FAMACHA score can be effective in improving parasite resistance because of its ease of recording, high heritability, and favorable genetic correlation with FEC. Additionally, runs of homozygosity islands related to parasite resistance could harbor important candidate genes for controlling this trait in Florida Cracker sheep.
Florida Cracker sheepgenetic parametersruns of homozygositygastrointestinal parasite resistanceheritage sheepsection-in-acceptanceAnimal Breeding and GeneticsIntroduction
Gastrointestinal parasites are a significant threat to sheep farming in the Southern United States; particularly, infection with Haemonchus contortus, which has detrimental effects on production, meat quality, and body condition scores in lambs (Zhong et al., 2016). One of the most sustainable strategies for controlling these infections is selecting for parasite-resistant animals. This method involves breeding individuals with a superior ability to regulate gastrointestinal parasites. Consequently, improved animal survival, health, and productive performance are expected due to cumulative and permanent genetic changes added to the new generations in the population. Additionally, resistant sheep secrete fewer H. contortus eggs in the feces, reducing pasture contamination (Woolaston and Baker, 1996).
Common phenotypes utilized for selecting parasite-resistant individuals include FAMACHA score (FAM), fecal egg count (FEC), packed cell volume (PCV), and productive traits (body weight, ADG, etc.) (Kaplan et al., 2004; Ngere et al., 2017; Ngere et al., 2018; Werne et al., 2023). Previous research on sheep populations suggests that the favored response to selection for parasite resistance originates from directional selection (Estrada-Reyes et al., 2019b), which creates an accumulation of runs of homozygosity (ROH) at specific loci (Chen et al., 2018; Islam et al., 2019). The ROH are contiguous segments of homozygous genotypes. These genotypes form part of identical by descent haplotypes transmitted by the parents to their offspring. The ROH are used to estimate inbreeding levels (Sams and Boyko, 2019) and are indicators of selection signatures (Kim et al., 2013; Metzger et al., 2015; Purfield et al., 2017). For this reason, ROH are commonly used to identify potential genomic regions under selection and candidate genes controlling economically important traits in livestock (Islam et al., 2019; He et al., 2020).
Florida Cracker sheep is a wool breed that evolved under natural selection and has good potential for meat production. The main features of this breed are its exceptional ability to control gastrointestinal parasites (Zajac et al., 1988; Amarante et al., 1999a; Amarante et al., 1999b; Estrada-Reyes et al., 2019a) and its adaptation to hot and humid conditions. This breed descends from Spanish sheep introduced to Florida in 1565 during the Spanish missions and the foundation of St. Augustine city (https://floridacrackersheep.com/history.html). Unfortunately, this breed is an endangered species under critical status (https://livestockconservancy.org/heritage-breeds/heritage-breeds-list/florida-cracker-sheep/). Hence, conservation and breeding programs are crucial to preserve this unique genetic resource. To rescue this breed, a combined breeding and conservation program is needed. This requires information about genetic parameters and inbreeding metrics; therefore, our objectives were 1) to estimate the heritability and genetic correlations for traits used to select for parasite resistance and 2) to study the distribution of ROH islands related to parasite resistance across the genome in the Florida Cracker sheep population.
Materials and methodsAnimal population and phenotype collection
The research protocol was approved by the University of Florida Institutional Animal Care and Use Committee (Approval number 201810108). The phenotypes used in this study were recorded on lambs between 3 and 5 months old during the summer from 2018 to 2020 at the largest commercial Florida Cracker sheep farm (Fairmeadow sheep farm) in Ocala, Florida. This farm uses a strict breeding program that avoids mating between related individuals. Before phenotype collection, experimental lambs were dewormed with levamisole (18 mg per kg of body weight). Then, ten days post-deworming, fecal samples were collected to verify FEC reduction. After this, lambs returned to the commercial farm and were naturally infected with H. contortus. The FEC assays were performed following the methodology described by Zajac and Conboy (2012). Briefly, fecal samples were collected from the rectum of each lamb and transported to the laboratory for analysis. Samples were processed using a saturated saline solution (NaCl). Quantification of eggs per gram of feces was performed by microscopic visualization using a McMaster chamber (Gordon and Whitlock, 1939).
Phenotypic records for FEC, packed cell volume (PCV), FAM, and body condition score (BCS, 1-5) were measured 38 days post-infection with H. contortus. The FEC assessment was performed as described above. The FAMACHA score was evaluated using the FAMACHA card (Van Wyk and Bath, 2002). The BCS was performed using a 1 − 5 scoring system and by palpation and visualization of the lumbar region according to the methodology described by Russel et al. (1969). The PCV was assessed using the microhematocrit method. PCV values from 27 − 45% were considered normal (Dargie and Allonby, 1975). The dataset contained 365 animals with phenotypic records at 38 days post-infection (natural Haemonchus contortus infection) for FEC, PCV, FAM, and BCS. The pedigree file comprised 695 animals born between 2016 and 2020, the progeny of 33 sires and 145 dams. The normality of data was assessed using the Shapiro-Wilk test. Then, the R “car” library was used to perform the Box-Cox transformation and to estimate the power parameter λ. Then, logarithmic transformation for FEC, which was not normally distributed, was carried out as follows: (log (FEC+1)).
Genotyping
Out of the 365 animals with phenotypic records, only 300 individuals with extreme phenotypes were selected for genotyping. Blood was collected from the jugular vein of each sheep using vacutainer tubes with anticoagulant EDTA to extract DNA using DNeasy Blood and Tissue Kit (Qiagen, Valencia, CA) according to the manufacturer’s instructions and stored at − 20°C. The DNA yield was calculated using the NanoDrop‐1000 at 260 nm (NanoDrop‐1000, Thermo Scientific). The DNA purity was evaluated using a ratio of 260/280 nm. Lambs were genotyped with the GGP Ovine 50K single nucleotide polymorphisms (SNPs) chip (GeneSeek, Inc., Lincoln, NE).
The initial SNP data included 45,205 SNPs and were utilized for quality control. Linkage Disequilibrium (LD) pruning was performed using the default options in the SNP and Variation Suite (SVS, Golden Helix, Inc., Bozeman, MT, USA), for this procedure, a sliding window of 50 SNPs with a 5 SNP increment, and a cutoff of R2 = 0.5 was utilized. Additional quality control was performed using the PREGSF90 software (Misztal et al., 2014) to remove duplicated genotypes, monomorphic markers, animals, and markers with call rate< 90%, markers with minor allele frequency< 0.05, markers with departure from the Hardy–Weinberg equilibrium (difference between expected and observed frequency of heterozygous) > 0.15, and animals with parent–progeny Mendelian conflicts. Additionally, non-mapped markers or markers mapped on sexual chromosomes were removed. After quality control, 38,429 SNP remained for 279 genotyped animals.
Genetic parameters estimation
Variance components, heritabilities, and genetic correlations for the FEC, PCV, FAM, and BCS were estimated using a multi-trait model with a Bayesian inference via the Gibbs sampling algorithm as implemented in the GIBBS3F90 program (Misztal et al., 2014). A single Gibbs chain with a total length of 300,000 iterations was generated. After discarding the initial 30,000 samples as burn-in, 1 in every ten samples was stored to compute means and standard deviations of the posterior distributions. The means were used as estimates of the genetic parameters, and their posterior standard deviations measured their estimation errors.
where FEC, PCV, FAM, and BCS are the traits as defined previously; y is the vector of phenotypic records; b is a vector of systematic effects: year (3 levels; 2018, 2019 and 2020), age group (3, 4 and 5 months old classes), and sex (2 levels; male and female); u is the vector of additive direct genetic effects; X, and Z are incidence matrices relating the elements of y to elements of b and u, respectively; e is the vector of random residual effects. All random effects were assumed to be multivariate, normally distributed with null expectations and the following covariance structure:
u~ N(0,Vu⊗H),ande~N(0,Re⊗I)
Vu is a 4 × 4 covariance matrix among traits for additive direct genetic effects, ⊗ is the Kronecker product, Re is a 4 × 4 residual error covariance matrix among traits, and H is a matrix that accounts for pedigree and genomic relationships among animals, and its inverse (H−1), which is needed to construct the mixed model equations, was obtained as in Aguilar et al. (2010).
H−1=A−1+[000G−1−A22−1]
where, A−1 is the inverse of the pedigree relationship matrix, A22−1 is the inverse of the pedigree relationship matrix among genotyped animals and, G−1 is the inverse of the genomic relationship matrix. The construction of the G matrix was based on VanRaden’s (2008) first method. Matrix G was blended with 5% of the block of the pedigree relationship matrix A corresponding to genotyped animals (A22) to avoid singularity problems. The rescaling of G to match A22 involved diagonals and off-diagonals as described in Chen et al. (2011).
Uniform prior distributions were assumed for b; b~p(b) where p(b)∝ constant. Scaled inverted Wishart prior distributions with scale matrix S and degrees of freedom parameter v, were assigned to covariance matrices for u and e to represent vague prior knowledge about these parameters, which presumably have little effect on results:
(Vu⊗H)~IW(S,v), S=[0] and v=−(p+1)=−5; p is equal to the dimension of Vu=4.
(Re⊗I)~IW(S,v), S and v are defined as above.
Heritability (h2) for each trait was computed using the following equation:
h2=σu2(σu2+σe2)
where σu2, and σe2 are the additive and residual variances, respectively.
Genetic correlations ( ρu) were calculated as:
ρu=σuT1,uT2(σuT12*σuT22)0.5
where σuT1,uT2 is the additive genetic covariance between the two traits, and σuT12 and σuT22 are the additive direct genetic variances for trait 1 and trait 2, respectively. Similarly, phenotypic correlations were computed.
Runs of homozygosity and inbreeding coefficient estimation
Runs of homozygosity and the inbreeding coefficient were estimated in 279 individuals and computed using the Golden Helix SNP and Variation Suite (SVS) 8.7.0 software (Golden Helix, Inc., Bozeman, MT, USA; www.goldenhelix.com). The following criteria were used to define the ROH: (1) one missing SNP and up to one possible heterozygous genotype were allowed in the ROH, (2) the minimum number of SNP that constituted a ROH was set to 30, (3) the minimum SNP density per ROH was set to one SNP every 100 Kb, and (4) the maximum gap between consecutive homozygous SNP was 250 Kb. The computed ROH were then classified into bins based on length:< 0.5 Mb, 0.5 – 2 Mb, > 2 – 5 Mb, > 5 – 10 Mb, > 10 – 20 Mb, and > 20 Mb.
The inbreeding coefficient was computed using the following formula:
FHOM=(O−E)/(L−E)
where FHOM is the probability of being homozygous by descent, O is the number of observed homozygotes for the individual over all markers, L is the number of genotyped autosomal markers, and E is the number of homozygous expected by chance. This inbreeding coefficient is equivalent to Wright’s within-subpopulation fixation index (Wright, 1922). Estimation of E was performed based on the following formula:
E=∑j=1L[1−2pjqjTAj(TAj−1)],
where TAj is the number of non-missing genotypes for marker j, and pj and qj are the allele frequencies p and q for marker j.
Incidence of common ROH and ROH island identification
The detection of genomic regions with high homozygosity was based on the number of times each SNP appeared in the ROH divided by the number of animals included in the analysis. These values were then plotted against the position of the SNP along the chromosome. ROH islands were analyzed using Golden Helix SVS, defined as clusters of ROH that were > 1000 Kb with a minimum of 30 SNPs and found in more than 20 samples.
Genomic coordinates of the ROH islands were utilized to identify the genes contained in these regions using the Ovis aries v4 in the Genome Data Viewer genome browser (https://www.ncbi.nlm.nih.gov/genome/?term=ovis+aries). The Quantitative Trait Locus (QTL) Database (https://www.animalgenome.org/cgi-bin/QTLdb/OA/index) was used to identify QTL regions that overlapped the ROH islands. The genes within each ROH island were further analyzed with the Panther Classification System (Mi et al., 2013) to identify gene ontology (GO) terms and significantly enriched genes with a False Discovery Rate< 0.05 using a ruminant model (Bos taurus) as reference. The Cytoscape v3.7.1 software was used for gene network graphics (Shannon et al., 2003).
ResultsHeritability estimates and genetic correlations
Descriptive statistics for all the phenotypes evaluated in this study are presented in Table 1. The mean for FEC was 1,427 (± 2,811) eggs per gram of feces and the mean for FAMACHA score was 2.75 (± 0.81). For PCV, the mean was 26.21(±6.73), and for BCS, the mean was 2.5 (±0.36). Table 2 presents the genetic and residual variances and heritability estimates for FEC, FAM, PCV, and BCS. Heritability estimates (± posterior standard deviation) were 0.33 ± 0.09 for FEC, 0.31 ± 0.10 for FAM, 0.22 ± 0.09 for PCV, and 0.19 ± 0.07 for BCS. The genetic correlations between the evaluated traits are presented in Table 3. The estimates for genetic correlations (± posterior standard deviation) among parasite resistance traits were 0.51 ± 0.21 (FEC-FAM), -0.74 ± 0.22 (FEC-PCV), and -0.18 ± 0.31 (PCV-FAM). Genetic correlations among parasite resistance traits and BCS were 0.32 ± 0.35 (FEC-BCS), -0.22 ± 0.31 (PCV-BCS), and -0.17 ± 0.28 (FAM-BCS); however, due to our small dataset, the posteriors standard deviations were high, yielding 95% central posterior intervals (Table 3) containing zero or with values out of the parameter space.
Descriptive statistics for fecal egg count (FEC), FAMACHA score (FAM), packed cell volume (PCV), and body condition score (BCS) in Florida Cracker sheep at 38 days post-infection naturally infected with H. contortus.
Trait
n
Mean (SD)
Min
Max
FEC (eggs/gram)
365
1,427 (2,811)
0
19,300
FAM (score)
365
2.75 (0.81)
1
5
PCV (%)
365
26.21 (6.73)
8.3
42
BCS (score)
365
2.5 (0.36)
1
3.25
Genetic (σ u2) and residual (σe2) variances, and heritability (h2) estimates for fecal egg count (FEC), FAMACHA score (FAM), packed cell volume (PCV), and body condition score (BCS) with 95% central posterior interval (in parentheses) in Florida Cracker sheep at 38 days post-infection naturally infected with H. contortus.
Trait
σ u2
σ e2
h2
FEC
0.39 (0.14,0.65)
0.78 (0.57,0.99)
0.33 (0.15,0.52)
FAM
0.20 (0.05,0.35)
0.44 (0.32,0.57)
0.31 (0.11,0.51)
PCV
5.39 (0.66,10.22)
18.35 (13.92,22.86)
0.22 (0.05,0.41
BCS
0.01 (0.00,0.03)
0.06 (0.05,0.07)
0.19 (0.05,0.32)
Genetic (above the diagonal) and phenotypic (below the diagonal) correlations among fecal egg count (FEC), FAMACHA score (FAM), packed cell volume (PCV), and body condition score (BCS) with 95% central posterior interval (in parentheses) in Florida Cracker sheep at 38 days post-infection naturally infected with H. contortus.
Trait
FEC
FAM
PCV
BCS
FEC
0.51 (0.09,0.94)
− 0.74 ( − 1.16, − 0.30)
0.32 ( − 0.37,1.02)
FAM
0.23 (0.13,0.34)
−0.18 ( − 0.79,0.44)
− 0.17 (− 0.73,0.40)
PCV
− 0.32 ( − 0.42, − 0.21)
− 0.33 ( − 0.44, − 0.23)
− 0.22 ( − 0.82,0.39)
BCS
− 0.16 (− 0.27, − 0.04)
− 0.25 ( − 0.36, − 0.14)
0.30 (0.20,0.40)
Runs of homozygosity and inbreeding coefficient
A total of 4,561 ROH ranging from 0.083 to 25.5 Mb were identified in the present study (Figure 1; Supplementary File 1). The average length of ROH was 3.33 Mb. The analysis of the distribution of ROH according to size (length) showed that the majority of the ROH were in the > 2 – 5 Mb category (47. 24%), followed by the 0.5 – 2 Mb category (26.17%). The longest ROH (> 20 Mb) only included 0.17% of the total ROH. The average inbreeding coefficient was − 0.02 ± 0.07.
Distribution of the number of runs of homozygosity (ROH) of different lengths (Mb) in Florida Cracker sheep.
Incidence of common ROH and ROH islands
The incidence of common ROH per SNP is presented in Figure 2 and Supplementary File 2. A higher incidence of ROH was observed in OAR5, OAR6, OAR12, OAR15, OAR19, and OAR21. The highest incidence of common runs per SNP was observed in OAR15, with over 120 incidences. The total number of ROH islands distributed across autosomal chromosomes was 113 (Figure 3). Detailed distribution of the ROH islands across the 26 ovine chromosomes is presented in Supplementary File 3. The number of islands ranged from a minimum of 1 cluster per chromosome (OAR16 and OAR23) to 9 clusters per chromosome (OAR2).
Incidence of common Runs of Homozygosity (ROH) per single nucleotide polymorphism (SNP) in Florida Cracker sheep.
Distribution of ROH islands per chromosome in Florida Cracker sheep.
Annotation of ROH islands
Gene annotation provided insight into the biological functions related to the ROH islands observed in Florida Cracker sheep. The genes identified in the ROH islands are presented in Supplementary File 4. The ROH islands in OAR1, 2, 11, 18, 19, and 20 included candidate genes for parasite resistance in Florida Cracker sheep and other sheep breeds in the US (IL12RB2, NFIL3, TRPM3, TNF, NOS2, CCR3, OVAR-DRA, and IL16) (Estrada-Reyes et al., 2018; Estrada-Reyes et al., 2019a; Estrada-Reyes et al., 2019b; Estrada-Reyes et al., 2021). The ROH islands overlapped parasite resistance and health QTLs previously reported in the Sheep QTL database (Supplementary File 4). The overlapped QTL included traits such as blood IgA level (OAR1, 4, 10, 12, 13, 15, 20, and 22), Haemonchus contortus FEC (OAR2, 3, 4, 10, 11, 12, 16, 20, 21 and 22), worm count (OAR2, 18, 22 and 26), strongyle FEC (OAR6, 7, 8, 9, 10, 14, 15, 17, 21, 22 and 26), Trichostrongylus adult and larva count (OAR2, 8 and 11), Trichostrongylus colubriformis FEC (OAR3, 11, 12 and 22), Nematodirus FEC (OAR14), hematocrit (OAR18, 22 and 26), and change in eosinophil number (OAR21). Most of the annotated genes in ROH islands were related to immune response processes (OAR2, 7, 8, 11, 12, 19, and 20) such as T-lymphocyte activation and proliferation, lymphocyte and leucocyte migration, granulocyte chemotaxis, activation of the innate immune response, immune cell extravasation, inflammatory response, chemokine, and cytokine signaling and production, and regulation of humoral immune response (Figure 4). The ROH islands identified in OAR7 included > 100 T-cell receptor genes, and OAR20 contained several (>14) major histocompatibility class I (MCH-I) and class II (MHC-II) genes (i.e., OVAR-DRB1, OVAR-DRB2, OVAR-DQB1, BOLA, etc.), and complement activation genes (C2, C4, C4A, CFB, BAG6, BTN1A1, NCR3, and BTNL2).
Gene network for enriched ROH islands controlling important immune response mechanisms in Florida Cracker sheep.
The PANTHER gene pathway analysis showed significantly enriched genes associated with the following biological processes: transcription factor binding (OAR1 and 3), kinase activity (OAR1), immune response processes (OAR2, 8, 11, 12, 19 and 20), homeostasis (OAR2), regulation of organelle organization (OAR2), DNA catabolism (OAR5), apoptosis (OAR5), cAMP signaling (OAR5), regulation of renin and endocrine system (OAR7), thrombin signaling (OAR7), regulation of cytosolic calcium and calcium-mediated signaling (OAR7 and OAR19), eye photoreceptor cell development (OAR7), hydro-lyase activity (OAR 9), linoleic acid metabolism and lipoxygenase activity (OAR11), sensory perception of smell (OAR 11 and 14), protein localization and modification (OAR12 and OAR26), regulation of nitrogen compound metabolism (OAR14 and 26), DNA damage signaling (OAR15 and 19), and nucleosome and chromosome organization (OAR20).
Discussion
Florida Cracker sheep have been naturally and artificially selected for resistance to gastrointestinal parasites (Estrada-Reyes et al., 2019a; Estrada-Reyes et al., 2021). This study estimated the genetic parameters for FEC, FAM, PCV, and BCS in Florida Cracker sheep naturally infected with H. contortus and the distribution of runs of homozygosity islands related to parasite resistance. The range and means for FAM and BCS recorded in our study were similar to those reported in other sheep populations, such as Merinos from South Africa (Dlamini et al., 2019; Snyman and Fisher, 2019). It is essential to highlight that Dlamini et al. (2019) studied a Merino population under selection for resistance to Haemonchus contortus since 2011; however, in this selected population, the mean FEC was higher (3,281 across several samplings) than in the Florida Cracker population (1,427). The h2 estimates for FEC and PCV (0.33 and 0.22, respectively) were of moderate magnitude; therefore, genetic progress could be feasible by selection for these traits in Florida Cracker sheep. However, the small population evaluated in this study might influence the accuracy of the genetic parameter estimates. Further monitoring of h2 estimates for parasite resistance is required for implementing breeding and conservation programs in this breed. Notter et al. (2018) reported heritability values of 0.19 and 0.24 for FEC in Katahdin lambs at weaning and post-weaning. Similar results were observed by Ngere et al. (2018) in Katahdin lambs, who reported heritability values for FEC ranging from 0.18 to 0.26 for weaning and from 0.23 to 0.46 for post-weaning. In this study, we evaluated Florida Cracker lambs during the post-weaning period, and our results are similar to the estimates reported by Ngere et al. (2018).
For the FAMACHA score, the heritability estimate (0.31) was higher when compared to reported values in Merino lambs [0.24; Riley and Van Wyk (2009)] and Dorper sheep [0.19; Ngere et al. (2017)]. However, lower h2 values for FEC (0.10), PCV (0.13), FAM (0.30), and BCS (0.17) have been observed in Santa Ines sheep, a Brazilian breed raised under hot and humid conditions (Rodrigues et al., 2021).
All the posterior standard deviations for the genetic correlations were high due to our study’s small number of records, highlighting the need for collecting more information. The 95% central posterior interval for all the correlations but the FEC-FAM correlation included zero or values out of the parameter space; thus, we could not make any conclusion about these estimates. The correlation FEC-FAM was 0.51 and the 95% central posterior interval (0.09,0.94). This correlation suggests that improvement in parasite resistance could be achieved in this population by selection based on FAM, a trait with additional advantages due to its high heritability and ease of recording; however, the ample 95% central posterior interval also remarks the need for more information to confirm this. Previous studies have reported higher estimates for the genetic correlations between FEC and FAMACHA scores ranging from 0.85 to 0.89 (Morris et al., 2004; Riley and Van Wyk, 2009; Mpetile et al., 2015) in sheep.
The FEC is considered the standard phenotype utilized for selection for parasite resistance (Woolastont and Piper, 1996; Morris et al., 1997; Woolastont and Windon, 2001). However, sheep producers do not commonly use FEC due to the associated costs. FAMACHA score has been successfully utilized as an alternative selection trait when H. contortus is the predominant gastrointestinal nematode (Kaplan et al., 2004) in sheep operations. This method only requires using a FAMACHA card to assess the level of anemia. Our results suggest that genetic selection for FAMACHA can be effective in improving FEC. Runs of homozygosity were used in this study to investigate potential genomic regions under selection. We also observed a low average inbreeding coefficient (− 0.02 ± 0.07). This negative estimate indicates that in this population there is, on average, more heterozygosity than expected, although the value does not differ significantly from zero. This inbreeding level is slightly lower than values reported for Sarda (0.04) and Sardinian Ancestral black (0.06) sheep from Italy by Cesarani et al. (2019). Similar estimates for average inbreeding were reported in Boer (0.01) and Nubian (0.01) goats (Hidalgo-Moreno et al., 2020). It is possible that the inbreeding coefficient reported in this manuscript could be a result of the breeding program utilized at the studied commercial farm. Analysis of ROH and subsequent PANTHER gene pathway analysis provided insight into the biological implications related to conserved regions found in the genome of the Florida Cracker sheep population. 113 ROH islands were distributed across the 26 ovine chromosomes in this breed. Most of these ROH islands were related to immune response processes such as T-lymphocyte activation and proliferation, lymphocyte and leucocyte migration, granulocyte chemotaxis, activation of innate immune response, immune cell extravasation, inflammatory response, chemokine, and cytokine signaling and production, and regulation of humoral immune response. Several of these mechanisms have been previously highlighted in Florida Cracker sheep infected naturally with H. contortus (Estrada-Reyes et al., 2019a; Estrada-Reyes et al., 2021). Previous studies with Florida Cracker sheep have focused on the identification of potential DNA variants (SNPs and CNVs) associated with parasite resistance (Estrada-Reyes et al., 2018; Estrada-Reyes et al., 2019a; Estrada-Reyes et al., 2021; Estrada-Reyes et al., 2022). The ROH islands identified in this manuscript overlapped these variants and other QTLs associated with FEC, adult worm count, IgA level, changes in eosinophil count and hematocrit, and included previously reported candidate genes with significant SNPs associated with FEC, red blood cell count, hemoglobin level, neutrophil count and PCV such as TRPM3, IL16 and TNF genes in this sheep breed (Estrada-Reyes et al., 2019a; Estrada-Reyes et al., 2021). Other important candidate genes included in these ROH islands were NOS2, CCR3, OVAR-DRA, and IL12RB2, which have been associated with parasite resistance in St. Croix and Katahdin sheep (Estrada-Reyes et al., 2018; Estrada-Reyes et al., 2019b; Estrada-Reyes et al., 2019c). Thus, it is possible that immunity-related genes within ROH islands and previously identified candidate genes, including TRPM3, IL16, TNF, IL12RB2, NOS2, and CCR3, may have a crucial role in controlling H. contortus and gastrointestinal parasites in this heritage breed.
The TRPM3 gene encodes a protein that forms a transient receptor potential channel. The transient receptor potential channels are calcium-permeable cation channels that allow Ca2+ influx and induce cell depolarization (Thiel et al., 2017). The exact role of the TRPM3 gene during H. contortus infections is still unknown, but this gene may play an essential role in controlling Ca2+ influx in immune cells during response processes against gastrointestinal parasites, including H. contortus. The IL16 gene regulates CD4+T lymphocyte production, chemotaxis, and macrophage polarization (Cruikshank et al., 2000; Huang et al., 2019). This was first proposed as a chemoattractant factor for T lymphocytes (Center and Cruikshank, 1982). The TNF is a cytokine that triggers inflammatory mechanisms, and it is commonly secreted by macrophages and neutrophils (Tamassia et al., 2019). Production of TNF by activated neutrophils (Zimmermann et al., 2015) and macrophages can mediate parasite resistance by inducing nitric oxide production during parasitic infections (Silva et al., 1995; Elmahallawy et al., 2021). In previous studies with Florida Cracker sheep, SNPs within TNF were associated with neutrophil count (Estrada-Reyes et al., 2019a). The IL12RB2 and NOS2 genes participate in signaling T helper 1 response (Thierfelder et al., 1996; Giordano et al., 2011). The CCR3 is present in eosinophils and is responsible for binding CCL11 or eotaxin, a specific eosinophil chemokine (Erin et al., 2002). This chemokine was also included in the ROH islands (OAR11) from our study. Further studies are required to investigate the role of these genes in Florida Cracker sheep production.
Conclusion
Parasite resistance traits (FEC, FAM, and PCV) are moderately heritable and suitable for selection in Florida Cracker sheep. The favorable genetic correlation between FAMACHA score and FEC suggest that genetic improvement for parasite resistance in this breed could be achieved by using FAMACHA score because of its ease of recording and reduced costs. Selection for parasite resistance results in the accumulation of ROH associated with this trait in Florida Cracker sheep. The ROH could contain important candidate genes for parasite resistance in this sheep breed. Further studies are required to validate our results.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Files, further inquiries can be directed to the corresponding author/s.
Ethics statement
The animal study was approved by University of Florida Institutional Animal Care and Use Committee (Approval number 201810108). The study was conducted in accordance with the local legislation and institutional requirements.
Author contributions
ZE-R designed the experiment. ZE-R, JH, IO, AP-C, TT, MI and GT assisted with data collection. JH and ZE-R conducted the data analysis. JH, ZE-R, IO, AP-C, TT, MI and GT drafted the manuscript. JH and ZE-R reviewed the final manuscript together. All authors contributed to the article and approved the submitted version.
Funding
The author(s) declare financial support was received for the research, authorship, and/or publication of this article. Neogene Corporation partially supported funding support. This study received partial funding from Neogen Corporation - GeneSeek Operations. The funder was not involved in the study design, collection, analysis, interpretation of data, the writing of this article or the decision to submit it for publication.
Acknowledgments
We acknowledge Carol Postley (Fairmeadow Farm), Dr. Owen Rae, and the students at the College of Veterinary Medicine at the University of Florida for their participation in phenotypic data collection. We gratefully acknowledge the valuable comments of Dr. Ignacy Misztal from the University of Georgia and the anonymous reviewers.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fanim.2023.1249470/full#supplementary-material