Advertisement

Genetic parameters and genomic regions associated with horn fly resistance in organic Holstein cattle

Published:September 02, 2021DOI:https://doi.org/10.3168/jds.2021-20366

      ABSTRACT

      Horn flies (Haematobia irritans [L.]) contribute to major economic losses of pastured cattle operations, particularly in organic herds because of limitations on control methods that can be used. The objectives of this research were to determine if resistance to horn flies is a heritable trait in organic Holstein cattle, determine associations with yield traits, and to detect genomic regions associated with fly infestation. Observations of fly load were recorded from 1,667 pastured Holstein cows, of which 640 were genotyped, on 13 organic dairies across the United States. Fly load score was determined using a 0 to 4 scale based on fly coverage from chine to loin on one side of the body, with 0 indicating few to no flies and 4 indicating high infestation. The scoring system was validated by counting flies from photographs taken at the time of scoring from 252 cows. To mitigate the effect of our data structure on potential selection bias effects on genetic parameter estimates, survival to subsequent lactations of scored animals and herd-mates that had been culled before the trial was accounted for as the trait stayability. Genetic parameters were estimated using single-step genomic analysis with 3-trait mixed models that included fly score, stayability, and a third phenotype. Model effects differed by variable, but fixed effects generally included a contemporary group, scorer, parity, and stage of lactation; random effects included animal, permanent environment, and residual error. A genome-wide association study was performed by decomposing estimated breeding values into marker effects to detect significant genomic regions associated with fly score. The rank correlation between the subjective fly score and the objective count was 0.79. The average heritability of fly score (± standard error) estimated across multiple models was 0.25 ± 0.04 when a known Holstein maternal grandsire was required and 0.19 ± 0.03 when only a known Holstein sire was required. Genetic correlation estimates with yield traits were moderately positive, but a greater fly load was associated with reduced yield after accounting for genetic merit. Lower fly loads were associated with white coat coloration; a significant genomic region on Bos taurus autosome 6 was identified that contains the gene KIT, which was the most plausible candidate gene for fly resistance because of its role in coat pattern and coloration. The magnitude of heritable variation in fly infestation is similar to other traits included in selection programs, suggesting that producers can select for resistance to horn flies.

      Key words

      INTRODUCTION

      The horn fly, Haematobia irritans [L.], is a member of the Muscidae family of flies that generally consists of parasites and pests to domestic animals and humans (
      • Moon R.D.
      Chapter 17: Muscid flies (Muscidae).
      ). As an obligate pest of cattle, the horn fly is considered the most prevalent cattle ectoparasite in North America. Annually, infestation results in over $700 million of lost profit in the United States with another $60 million spent by producers on fly control (
      • Byford R.L.
      • Craig M.E.
      • Crosby B.L.
      A review of ectoparasites and their effect on cattle production.
      ;
      • Foil L.D.
      • Hogsette J.A.
      Biology and control of tabanids, stable flies and horn flies.
      ).
      Organic cattle producers face some of the greatest challenges associated with horn flies. In a recent survey of organic dairy producers in Ohio, about half of the respondents reported that fly control was a major management challenge in their herds (
      • Brock C.C.
      • Pempek J.A.
      • Jackson-Smith D.
      • Weaver K.
      • da Costa L.
      • Habing G.G.
      Organic dairy producer experiences and decisions related to disease prevention and treatment.
      ). In the United States, organic livestock are required to have access to the outdoors, when weather permits (). As horn flies are most prolific in pasture settings, organic operations are particularly susceptible to having infested cattle (
      • Moon R.D.
      Chapter 17: Muscid flies (Muscidae).
      ). Further, organic herds in the United States must abide by the National List of Allowed and Prohibited Substances (), which prohibits the use of all synthetic chemical pest control in livestock. Organic regulations also forbid the use of antibiotics in livestock, which means when an organic animal contracts a bacterial infection spread by horn flies, it cannot be treated with traditional drugs. Thus, the prevention of horn fly infestation in organic herds is essential.
      A proposed method to accommodate these needs is breeding cattle to be naturally resistant to horn flies (
      • Brown Jr., A.H.
      • Steelman C.D.
      • Johnson Z.B.
      • Rosenkrans Jr., C.F.
      • Brasuell T.M.
      Estimates of repeatability and heritability of horn fly resistance in beef cattle.
      ;
      • Kerlin R.L.
      • Allingham P.G.
      Acquired immune response of cattle exposed to buffalo fly (Haematobia irritans exigua).
      ). Cattle that consistently host fewer horn flies than their herd-mates are considered to have resistance to horn flies. Genetic parameters of fly resistance have been estimated in beef breeds and beef producers have indicated that they would pay a premium for bulls with genetic resistance to horn flies (
      • Brown Jr., A.H.
      • Steelman C.D.
      • Johnson Z.B.
      • Rosenkrans Jr., C.F.
      • Brasuell T.M.
      Estimates of repeatability and heritability of horn fly resistance in beef cattle.
      ;
      • Fraga A.B.
      • Mello De Alencar M.
      • Andrade De Figueiredo L.
      • Razook A.G.
      • Noely J.
      • Cyrillo S.G.
      Genetic analysis of the infestation of females of the Caracu cattle breed by Horn Fly (Haematobia irritans irritans) (L.) (Diptera, Muscidae).
      ;
      • McKay L.
      • DeLong K.L.
      • Schexnayder S.
      • Griffith A.P.
      • Taylor D.B.
      • Olafson P.
      • Trout Fryxell R.T.
      Cow-calf producers' willingness to pay for bulls resistant to horn flies (Diptera: Muscidae).
      ). However, this trait has not been investigated extensively in dairy breeds.
      Therefore, the objectives of this study were to estimate the genetic parameters of horn fly resistance, determine traits correlated with horn fly resistance, and to determine candidate genes associated with horn fly resistance in organic Holstein cows.

      MATERIALS AND METHODS

      Animals

      All protocols were approved by the Pennsylvania State University Institutional Animal Care and Use Committee (ID: PRAMS201747560). Dairy cattle were located across the contiguous United States on commercial organic herds (n = 11) and organic university research herds (n = 2). Herds had varying levels and methods of fly control implemented including USDA certified organic sprays and drenches and walk-through fly vacuums. Fly load scores were visually assessed and recorded by trained scorers on lactating cows at least once annually over the years 2017 to 2019. Tail hair was collected from a random subset of animals to facilitate genomic analysis; cows were genotyped for either 19K (Zoetis Clarifide) or 42K (Neogen GGP Dairy) markers. Genotypes were submitted to the Council on Dairy Cattle Breeding and imputed to the 79K marker genotype used in national genetic evaluations (
      • VanRaden P.M.
      • O'Connell J.R.
      • Wiggans G.R.
      • Weigel K.A.
      Genomic evaluations with many more genotypes.
      ;
      • Wiggans G.R.
      • VanRaden P.M.
      • Null D.J.
      • Cole J.B.
      Genomic predictions using more markers and gene tests.
      ).

      Scoring System

      To identify the degree of horn fly resistance a cow had, a method of quantification of fly infestation was developed. A scoring system similar to that of
      • Fraga A.B.
      • Mello De Alencar M.
      • Andrade De Figueiredo L.
      • Razook A.G.
      • Noely J.
      • Cyrillo S.G.
      Genetic analysis of the infestation of females of the Caracu cattle breed by Horn Fly (Haematobia irritans irritans) (L.) (Diptera, Muscidae).
      was used to quantify fly load. Fly density was visually assessed based on fly coverage from chine to loin on one side of the body. Fly scores were assigned on a scale of 0 to 4; a score of 0 indicated few to no flies, whereas a score of 4 indicated high infestation. The scoring system was validated from photographs taken at the time of scoring on a subset of observations. Flies in the photographs located in the same region scored were counted digitally using the multipoint feature of ImageJ (
      • Schneider C.A.
      • Rasband W.S.
      • Eliceiri K.W.
      NIH image to ImageJ: 25 years of image analysis.
      ; Figure 1). Fly scores were assigned to grazing animals between the months of May and September, when horn fly prevalence is greatest (
      • Lysyk T.J.
      Effect of larval rearing temperature and maternal photoperiod on diapause in the horn fly (Diptera: Muscidae).
      ). When visiting a farm to assign fly scores, our objective was to efficiently score many cows rather than scoring any one particular cow. One or more trained evaluators moved around the pasture or barn and recorded cow identification number and scored all cows within the vicinity; personnel selected different areas of the pasture to work when multiple evaluators were present at the same farm. Some cows may have been scored by multiple evaluators, could be scored multiple times by the same evaluator, or could have evaded scoring in some instances due to movement in the pasture.
      Figure thumbnail gr1
      Figure 1Photographs of individuals that were assigned fly scores of 0 through 4, respectively. The yellow represents where a horn fly was counted. Images are labeled with fly score assigned in the top left corner and total fly count in the bottom left corner.

      Evaluating Coat Color

      Cattle color composition was determined from the photographs that were taken at the time of fly scoring. Color composition was estimated by trained evaluators who assigned percentage values of the colors black, white, and red based on the photograph. Color evaluators were blinded to the fly load components of the study. If an individual animal was photographed multiple times, the color composition estimates were averaged across photos. Most of the Holsteins in our data were black and white so observations of red percentage were not included in analyses. This resulted in white coat color composition values being the inverse values of black coat color composition. As such, all analyses of coat color resulted in nearly identical and, regarding correlations, inverted parameter estimates. Thus, only analyses of white coat color composition are reported.

      Stayability

      Stayability was considered the likelihood of an animal staying within a herd up to a certain time point (
      • Hudson G. F.S.
      • Van Vleck L.D.
      Relationship between production and stayability in Holstein cattle.
      ). Death and culling, both voluntary and involuntary, were considered factors for leaving the herd. Management on the organic herds involved in this study varied immensely; not all herds were on DHIA testing nor exclusively used commercially available AI sires. This led to incomplete records within herds and record inconsistency between herds, which could have biased genetic parameter estimates if fly score was correlated with survival earlier in life (
      • Árnason T.
      • Albertsdóttir E.
      • Fikse W.F.
      • Eriksson S.
      • Sigurdsson Á.
      Estimation of genetic parameters and response to selection for a continuous trait subject to culling before testing.
      ). Further, the use of natural service bulls without genetic evaluations in some herds restricted pedigree depth for animals descending from such bulls. In our data, cows in their second or later lactation at the start of the study did not have records from first lactation and may have had contemporaries that left the herd before scoring. For these reasons stayability was used in multivariate analyses to mitigate any effects of selection bias on our genetic parameter estimates (
      • Hardie L.C.
      • Heins B.J.
      • Dechow C.D.
      Genetic parameters for stayability of Holsteins in US organic herds.
      ).
      Historic herd records maintained by the producers were used to compile a comprehensive pedigree and to determine stayability. Stayability was treated as a binary variable that was quantified by parity (1 = remained in the herd through the completion a given parity; 0 = otherwise). For example, if a cow that was scored during her first and second lactation was culled before 305 d in her second lactation, she would receive a stayability score of 1 for her first parity and 0 for her second. Cows that had yet to be given the opportunity to stay in the herd through fifth lactation were assigned missing stayability values for their future parities.

      Data Management

      Any records beyond fifth lactation were removed from the data. If phenotypes were evaluated on an individual multiple times by the same scorer on the same day, numeric scores were averaged for analysis. Scores assigned to the same cow by multiple evaluators on the same day were retained from only the evaluator that had scored the most animals through the duration of the study. Best predictions (
      • Cole J.B.
      • VanRaden P.M.
      A Manual for Use of BESTPRED: A Program for Estimation of Lactation Yield and Persistency Using Best Prediction.
      ) of the mature equivalent of a 305-d lactation (ME305) were generated from available monthly DHIA test records for milk, fat, and protein yield, and SCS, which were merged with the first instance of scoring observation within a given parity. As not all herds used DHIA testing, ME305 information was not available for all cows that received fly scores or stayability records.
      Breed composition of individual cows was calculated from their pedigrees; however, genomic breed base representation (BBR) was used for individuals that were genotyped (
      • VanRaden P.M.
      • Cooper T.A.
      Genomic evaluations and breed composition for crossbred U.S. dairy cattle.
      ;
      • Norman H.D.
      • VanRaden P.M.
      • Megonigal J.H.
      • Dürr J.W.
      • Cooper T.A.
      Breed base representation in dairy animals of five breeds.
      ). Initial records had diverse breed representation but still primarily consisted of Holsteins. To restrict our data to Holstein lineage, genotyped animals were excluded if their BBR was less than 87.5%. For inclusion of animals with only pedigree information, Holstein-sired animals were included that had no other breeds identified in their 4-generation pedigree. Data were defined as 2 sets: animals with 75% or greater known Holstein lineage (HO75) and those with a known Holstein sire (HO50) and no known ancestors of other breeds. The HO75 data set was the primary focus for genetic analyses with HO50 used to facilitate greater statistical power from models that included fly count and coat color, as those traits had relatively low sample sizes. Table 1 summarizes the available phenotype records for HO75 and HO50 used for analysis. The HO75 data set contained 1,639 records of fly score from 1,146 individuals that were scored at least once during their first to fifth lactation and 4,743 records of ME305 traits from 2,398 individuals from lactations 1 to 5. A total of 52,351 records of stayability for lactations 1 to 5 were available from 13,529 individuals in HO75. In the HO50 data, the same lactation criteria resulted in 2,231 available records of fly score from 1,667 individuals and 140,445 records of stayability from 35,395 individuals. A subset of the assigned fly scores in the HO75 data also had counted flies (n = 281) from pictures of 229 individual animals; in the HO50 data 305 fly counts were available from 252 individuals. In HO75 and HO50 there were 390 and 503 observations of coat color composition, respectively. As counting flies was an inefficient process, horn flies were only counted from enough photographs to validate our scoring system, whereas coat color composition was evaluated from all available photographs. Thus, there were more observations of coat color composition than fly count. Genotypes were available from 2,273 individuals that had minimum Holstein BBR of 87.5%, of which 640 had a fly score and 234 had coat color evaluated.
      Table 1Descriptive statistics of traits considered for genomic analysis
      Data set
      HO75 = animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree; HO50 = animals with known Holstein sire and no other known breeds in 4-generation pedigree.
      Trait
      ME305 = mature equivalent of a 305-d lactation.
      Number of observationsNumber of cowsMean ± SD
      HO75Fly score1,6391,1461.02 ± 0.92
      Fly count (raw)28122944 ± 69
      Fly count (ln)2812293.06 ± 1.15
      Stayability52,35113,5290.35 ± 0.48
      White coat color (%)39039036.7 ± 32.0
      ME305 milk yield (kg)4,7432,3989,572 ± 1,710
      ME305 fat yield (kg)4,7432,398373 ± 67.4
      ME305 protein yield (kg)4,7432,398309 ± 50.6
      ME305 SCS4,7432,3982.85 ± 1.13
      HO50Fly score2,2311,6670.99 ± 0.91
      Fly count (raw)30525242 ± 67
      Fly count (ln)3052523.02 ± 1.14
      Stayability140,44535,3950.36 ± 0.48
      White coat color (%)50350335.4 ± 31.4
      1 HO75 = animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree; HO50 = animals with known Holstein sire and no other known breeds in 4-generation pedigree.
      2 ME305 = mature equivalent of a 305-d lactation.

      Comparison of Fly Score and Count

      To validate the scoring system, the Spearman rank correlation between the flies counted from photographs (n = 369) and fly score assigned at the time of photographing was determined using the CORR procedure of SAS (v. 9.4 SAS Institute Inc.) for the HO50 data. For this analysis, fly scores and fly counts from photographs acquired on the same date were not averaged, as we were interested in the correlation of the score with the unique number of flies counted at the time of scoring. The median and range of the number of flies for each fly score category were determined. For genetic analyses, fly count was transformed to a natural log value because the raw counts were not normally distributed.

      Genomic Analysis

      A single-step genomic analysis was performed using the methods of
      • Aguilar I.
      • Tsuruta S.
      • Legarra A.
      PREGSF90-POSTGSF90: Computational tools for the implementation of single-step genomic selection and genome-wide association with ungenotyped individuals in BLUPF90 programs. Proc. 10th World Congr. Genet. Appl. Livest. Prod. American Society of Animal Science.
      and implemented in the BLUPF90 family of programs (
      • Misztal I.
      • Tsuruta S.
      • Lourenco D.
      • Aguilar I.
      • Legarra A.
      • Vitezica Z.
      Manual for BLUPF90 family of programs.
      ). Models were first analyzed considering only pedigree relationship with the program ASREML (
      • Gilmour A.R.
      • Gogel B.J.
      • Cullis B.R.
      • Welham S.J.
      • Thompson R.
      ASReml User Guide Release 4.1.
      ). The genetic parameters estimated by the pedigree analysis were then used as starting values in our single-step genomic analysis. Genetic and phenotypic correlations were generally in the same direction as those estimated by genomic analysis with reduced standard error and fly score heritability estimates were slightly increased with smaller standard errors. As such, only results from single-step genomic analysis are reported. The inverse of the pedigree relationship matrix (A−1) was blended with the inverse of the genomic relationship matrix (G−1) of genotyped animals into a single inverse relationship matrix (H−1). A series of the following 3-trait models were used for multivariate genetic analysis:
      yi=Xiβi+Zaiai+Zpipi+Zh2h2+ei,


      where y is a vector observations for trait i which includes fly score (i = 1), stayability (i = 2), and one additional trait (i = 3); X is an incidence matrix of fixed effects; βi is a vector of fixed effects corresponding to trait i; Zai is the incidence matrix of animal effects and ai is the vector of random animal effects for trait i; Zpi and pi are the incidence matrix and random permanent environment effect of trait i, respectively; Zh2 and h2 are the incidence matrix and random herd-year-season of birth by parity effects for trait 2, stayability; and ei = the random effect of trait i.
      All models included fly score, stayability, and one of the following: natural log of fly count, ME305 of milk yield (kg), ME305 fat yield (kg), ME305 protein yield (kg), ME305 SCS, or percentage of white coat coloration. Tertiary traits were primarily included in our models to determine the genetic correlation between fly resistance and production and health traits. Fixed and random effects varied by trait. For fly score and fly count fixed effects included scorer, the interaction of date scored and herd, stage of lactation by week, and parity. Herd-year-season of parity and herd-year-season of birth were both considered for the model as well, but subclasses of those traits were small and did not significantly contribute to fly score variation. Instead, the fixed effect of date scored by herd accounted for herd-to-herd variation, variation caused by seasonality of scoring, and temporal effects of any fly control methods that were used. The model also included a random animal effect. Despite the substantial repeated measures of fly score, models run with fly score or fly count would not converge when a random permanent environment effect was included. As such, the permanent environment effect was assumed to be 0 from our data and the effect was excluded from fly score and fly count models.
      For stayability, fixed effects were parity and herd-year-season of birth. The interaction of the 2 fixed effects, herd-year-season of birth and parity, was included as a random effect because in some instances their subclasses had very few observations. Permanent environment and animal were also random effects in the stayability model. For ME305 traits fixed effects included parity, herd-year-season of calving, and age at calving, and permanent environment and animal were random effects. As percentage of white coat coloration remains approximately unchanged over lactations, a simpler model with the fixed effect of herd-year-season of birth and random effect of animal was used for color composition. Each model was run on the HO75 data, and all models but those including ME305 traits were run on the HO50 data set. Convergence was not achieved in either the HO75 nor the HO50 3-trait models of fly score, stayability, and fly count. Thus, the same model was fit using only pedigree-relatedness in ASREML (
      • Gilmour A.R.
      • Gogel B.J.
      • Cullis B.R.
      • Welham S.J.
      • Thompson R.
      ASReml User Guide Release 4.1.
      ).
      The AIREMLF90 program was used to estimate the genetic variance of each trait; estimated parameters were then used by BLUPF90 to generate EBV and SNP solutions (
      • Misztal I.
      • Tsuruta S.
      • Lourenco D.
      • Aguilar I.
      • Legarra A.
      • Vitezica Z.
      Manual for BLUPF90 family of programs.
      ). Out of the 79,060 genetic markers, 73,560 were included in the analysis after quality control measures removed monomorphic SNP (n = 45), SNP with call rates below 0.90 (n = 3,011), and SNP with minor allele frequencies below 0.05 (n = 2,588). Further, 47 animal genotypes were removed due to parent-progeny Mendelian conflicts, resulting in a total of 2,226 genotyped animals in the data.
      A GWAS was then performed to decompose the genomic EBV of fly score and percentage of white coloration into SNP effects by the program POSTGSF90 (
      • Misztal I.
      • Tsuruta S.
      • Lourenco D.
      • Aguilar I.
      • Legarra A.
      • Vitezica Z.
      Manual for BLUPF90 family of programs.
      ). The 3-trait model of HO75 fly score, stayability, and percentage of white coloration was used for the GWAS presented based on the phenotypic correlation estimate between fly score and white coloration. Significance levels for each SNP and variance explained by 1-Mb regions were determined. A Bonferroni correction was used to account for repeated testing; for a desired significance level of 0.05, a SNP was considered significant at P < 0.00000068. Genes encompassing the most significant 1-Mb region and those that contained the most significant SNP were identified using the Bos taurus ARS-UCD1.2 reference genome available via the National Center for Biotechnology Information and their search tool, GenBank (
      • O'Leary N.A.
      • Wright M.W.
      • Brister J.R.
      • Ciufo S.
      • Haddad D.
      • McVeigh R.
      • Rajput B.
      • Robbertse B.
      • Smith-White B.
      • Ako-Adjei D.
      • Astashyn A.
      • Badretdin A.
      • Bao Y.
      • Blinkova O.
      • Brover V.
      • Chetvernin V.
      • Choi J.
      • Cox E.
      • Ermolaeva O.
      • Farrell C.M.
      • Goldfarb T.
      • Gupta T.
      • Haft D.
      • Hatcher E.
      • Hlavina W.
      • Joardar V.S.
      • Kodali V.K.
      • Li W.
      • Maglott D.
      • Masterson P.
      • McGarvey K.M.
      • Murphy M.R.
      • O'Neill K.
      • Pujar S.
      • Rangwala S.H.
      • Rausch D.
      • Riddick L.D.
      • Schoch C.
      • Shkeda A.
      • Storz S.S.
      • Sun H.
      • Thibaud-Nissen F.
      • Tolstoy I.
      • Tully R.E.
      • Vatsan A.R.
      • Wallin C.
      • Webb D.
      • Wu W.
      • Landrum M.J.
      • Kimchi A.
      • Tatusova T.
      • DiCuccio M.
      • Kitts P.
      • Murphy T.D.
      • Pruitt K.D.
      Reference sequence (RefSeq) database at NCBI: Current status, taxonomic expansion, and functional annotation.
      ;
      • Clark K.
      • Karsch-Mizrachi I.
      • Lipman D.J.
      • Ostell J.
      • Sayers E.W.
      GenBank.
      ).
      The genetic parameters of heritability, repeatability, genetic correlation, and phenotypic correlation were calculated using the methods of
      • Falconer D.S.
      • Mackay T. F.C.
      Introduction to Quantitative Genetics.
      . The equation to estimate heritability was h2=σA2/σP2, where h2 is the narrow sense heritability of the respective trait; σa2 is the genetic variance for the respective trait; and σP2 is the phenotypic variance which is the sum of genetic variance, permanent environmental variance, and residual variance for the respective trait.
      Repeatability was estimated as r=σPA2/σP2, where r is the repeatability of the respective trait; σPA2 is the producing ability, or an individual's performance potential for the respective trait, which is the sum of genetic variance and permanent environmental variance; and σP2 is the phenotypic variance. Because repeatability was calculated this way, traits with permanent environmental variance of 0 had repeatability equivalent to heritability.
      Genetic correlation was estimated as rg=σAXAy/σAx2σAy2, where rg represents the genetic correlation between 2 traits; σAXAy is the genetic covariance of the 2 traits; and σAx2 and σAy2 are the genetic variances of the respective traits. Phenotypic correlation was estimated by rp=σPXPy/σPx2σPy2, where rp represents the phenotypic correlation between 2 traits; σPXPy is the phenotypic covariance of the 2 traits which is the sum of genetic covariance, permanent environmental covariance, and residual covariance of the 2 traits; and σPx2 and σPy2 are the phenotypic variances of the respective traits. The standard errors of each genetic parameter were calculated using the methods of
      • Houle D.
      • Meyer K.
      Estimating sampling error of evolutionary statistics based on genetic covariance matrices using maximum likelihood.
      through Monte Carlo sampling of distributions (n = 10,000) generated from the estimated variance components of each trait and the sampling variance matrix. The sampling variances of genetic parameter estimates of fly score and stayability across models were minimal (<0.001 for each parameter). Parameter estimates from each multivariate model were averaged across analyses for the respective data sets to best represent the multiple analyses performed.

      Fly Score and Test-Day Yield

      Because ME305 data reflects total lactation yield including seasons with minimal horn fly infestation, the effect of fly load on short-term production may be unclear. Therefore, we merged fly scores with test-day records (n = 1,197) of 886 animals that occurred within 2 wk of fly assessment in the HO75 data set. There were 1,176 records of milk yield, and fly score from 876 cows, and 954 records of fat and protein yield from 661 cows. The EBV from the genomic analysis of yield described above was also included in the data set to determine effects of fly load that are independent of the cow's genetic merit for yield.
      Test-day records were evaluated using the MIXED procedure of SAS with the following full model:
      yijklmn=μ+Flyi+be×EBV+Lactj+Datek+l=15bl×DIMl×LGm+Cown+Datek×LGm+εijklmn,


      where y = test-day milk yield, fat yield, or protein yield; Fly = fly score i (1 to 3 where scores of 3 and 4 were collapsed into a single category); µ = the intercept; be = coefficient of regression on EBV for trait i; Lact = lactation j (1 to 5); Date = date scored k; bl = coefficient of regression on the interaction of DIM through the fifth order polynomial of DIM and lactation group m (LG; 1 to 3 where lactations 3, 4, and 5 were collapsed into a single category); Cow = the random effect of individual cow n; the random effect of the interaction of date k with lactation group m; and ε = random error. Least squares means (LSM) were estimated by fly score with significance declared at P < 0.05, and trends declared at P < 0.10. Reduced models that excluded the effect of the respective EBV were also evaluated.
      In the fat yield model, the interaction of DIM and lactation group was only significant through the first order polynomial of DIM, so all others were eliminated. Similarly, in the test-day protein yield model DIM interaction was eliminated following the third-order polynomial level.

      Color Analyses

      The GWAS identified Hapmap27516-BTC-042465 as one of the SNP that significantly influenced fly score, which was located in a candidate gene, KIT proto-oncogene, receptor tyrosine kinase (KIT). Previously established relationships between KIT and piebaldism (
      • Grosz M.D.
      • MacNeil M.D.
      The “spotted” locus maps to bovine chromosome 6 in a Hereford-cross population.
      ;
      • Reinsch N.
      • Thomsen H.
      • Xu N.
      • Brink M.
      • Looft C.
      • Kalm E.
      • Brockmann G.A.
      • Grupe S.
      • Kühn C.
      • Schwerin M.
      • Leyhe B.
      • Hiendleder S.
      • Erhardt G.
      • Medjugorac I.
      • Russ I.
      • Fö M.
      • Reents R.
      • Averdunk G.
      A QTL for the degree of spotting in cattle shows synteny with the KIT locus on chromosome 6.
      ) led us to investigate the impact of variants of Hapmap27516-BTC-042465 and percentage of white coat coloration on fly score.
      Genotyped Holsteins with BBR greater than or equal to 87.5% and available color observations (n = 234) were included in the following analyses. To determine the relationship between the coat color of Holstein cows and variants of Hapmap27516-BTC-042465, proportion of white coloration was evaluated with the following linear model using the MIXED procedure of SAS: yi=μ+SNPi+εi, where y = proportion of white coat coloration; SNP = Hapmap27516-BTC-042465 variant i (0, 1, or 2); and ε = random error. Least squares means were estimated for SNP variant effect with significance declared at P < 0.05, and trends declared at P < 0.10.
      The same animals were then merged with available fly scores, reducing the total animals to 233 but increasing phenotypes as there were 435 fly scores among them. The MIXED procedure of SAS was used to evaluate the association between fly score and variants of Hapmap27516-BTC-042465 in the following model:
      yijklm=μ+SNPi+Herdj×Datek+Scorerl+Cowm+εijklm,


      where y = fly score; SNP = Hapmap27516-BTC-042465 variant i (0, 1, or 2); Herd = herd animal was scored in j; Date = date scored k; Scorer = individual that assigned fly score l; Cow = the random effect of individual cow m; and ε = random error. The same model was run again with the addition of percentage of white coloration as a fixed effect. For all models, LSM were estimated for the effect of SNP variant.

      RESULTS

      Fly Score Verification

      Fly Score was significantly correlated with fly count (ρ = 0.79; P < 0.0001). Cows had a range of 0 to 600 horn flies counted on a single side of their body; the median number of horn flies was 16. Table 2 displays descriptive statistics of fly count by fly score. Cows that had a fly score of 0 had a mean of 6 flies and median of 5 flies counted on one side of the body; cows that were scored 1, 2, 3, or 4, had mean fly counts of 25, 67, 145, and 351 and median fly counts of 19, 59, 127, and 318, respectively. A pattern of increasing mean and median fly counts was observed. Figure 2 displays the linear relationship between the 2 variables.
      Table 2Descriptive statistics of fly counts by fly score
      Fly scoreNo. of observationsMeanMedianRange
      0114650–33
      117625190–192
      252675910–283
      32114512728–354
      46351318194–600
      Figure thumbnail gr2
      Figure 2Fly score by fly count. Fly scores assigned in the field are represented on the x-axis and the fly count from photographs of individuals at scoring are in the y-axis.

      Genetic Parameter Estimates

      Heritability and repeatability estimates for evaluated traits are presented in Table 3. The values for fly scores were averaged across models resulting in an estimated fly score heritability of 0.25 ± 0.04 for HO75. The heritability estimate of fly score (0.19 ± 0.03) from the HO50 data set was somewhat less than for HO75, possibly reflecting more accurate and complete pedigrees in HO75. Models of fly count did not converge in the genomic analysis. Averaged stayability heritability was estimated as 0.06 ± 0.01 and repeatability was estimated as 0.71 ± 0.004 from HO75, respective estimates were 0.09 ± 0.01 and 0.71 ± 0.002 from HO50.
      Table 3Heritability and repeatability estimates of traits from single-step genomic analysis
      Data set
      HO75 = animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree; HO50 = animals with known Holstein sire and no other known breeds in 4-generation pedigree.
      TraitHeritability ± SERepeatability ± SE
      HO75Fly score0.25 ± 0.04
      Repeatability was not calculated because permanent environment effects were assumed to be 0.
      Stayability0.06 ± 0.010.71 ± 0.004
      Fly countNC
      NC = model did not converge.
      NC
      ME305 milk yield0.15 ± 0.030.43 ± 0.02
      ME305 fat yield0.17 ± 0.030.38 ± 0.02
      ME305 protein yield0.12 ± 0.030.47 ± 0.02
      ME305 SCS0.09 ± 0.030.42 ± 0.02
      Proportion of white coloration0.46 ± 0.20
      HO50Fly score0.19 ± 0.03
      Fly countNCNC
      Stayability0.09 ± 0.010.71 ± 0.002
      Proportion of white coloration0.44 ± 0.18
      1 HO75 = animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree; HO50 = animals with known Holstein sire and no other known breeds in 4-generation pedigree.
      2 Repeatability was not calculated because permanent environment effects were assumed to be 0.
      3 NC = model did not converge.
      Genetic and phenotypic correlations estimated by multitrait models are presented in Table 4. The averaged genetic correlations between fly score and stayability from HO50 and HO75 were moderate (0.16 ± 0.14 and 0.25 ± 0.17, respectively). A small, positive phenotypic correlation was estimated between the 2 traits from the averaged HO50 models and averaged HO75 models (0.05 ± 0.03 and 0.09 ± 0.03, respectively). Phenotypic correlations between fly load and coat coloration were significant in both data sets. Proportion of white coat coloration was negatively genetically (−0.81 ± 0.60 and −0.73 ± 0.78 in HO50 and HO75, respectively) and phenotypically (−0.34 ± 0.08 and −0.43 ± 0.09 in HO50 and HO75, respectively) correlated with fly score. Significant positive genetic correlations were estimated between fly load and ME305 milk (0.35 ± 0.14), fat (0.29 ± 0.13), and protein (0.44 ± 0.16). The ME305 SCS and fly load were not significantly correlated.
      Table 4Genetic and phenotypic correlations ± SE of fly score and stayability with production and other traits based on single-step genomic analysis
      HO75 = animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree; HO50 = animals with known Holstein sire and no other known breeds in 4-generation pedigree; ME305 = mature equivalent of a 305-d lactation.
      Data setTraitCorrelation typeFly scoreStayability
      HO75Fly scoreGenetic0.25 ± 0.17
      StayabilityPhenotypic0.09 ± 0.03
      SE is less than half of the correlation.
      Fly countGenetic NC
      NC = model did not converge.
      NC
      PhenotypicNCNC
      ME305 milk yieldGenetic0.35 ± 0.14
      SE is less than half of the correlation.
      0.50 ± 0.20
      SE is less than half of the correlation.
      Phenotypic0.07 ± 0.04−0.10 ± 0.02
      SE is less than half of the correlation.
      ME305 fat yieldGenetic0.29 ± 0.13
      SE is less than half of the correlation.
      0.17 ± 0.18
      Phenotypic0.06 ± 0.04−0.13 ± 0.02
      SE is less than half of the correlation.
      ME305 protein yieldGenetic0.44 ± 0.16
      SE is less than half of the correlation.
      0.58 ± 0.22
      SE is less than half of the correlation.
      Phenotypic0.06 ± 0.04−0.08 ± 0.02
      SE is less than half of the correlation.
      ME305 SCSGenetic0.08 ± 0.20−0.34 ± 0.24
      Phenotypic−0.004 ± 0.04−0.10 ± 0.02
      SE is less than half of the correlation.
      White coloration %Genetic−0.73 ± 0.78−0.23 ± 0.40
      Phenotypic−0.43 ± 0.09
      SE is less than half of the correlation.
      −0.08 ± 0.06
      HO50Fly scoreGenetic0.16 ± 0.14
      StayabilityPhenotypic0.05 ± 0.03
      Fly countGeneticNCNC
      PhenotypicNCNC
      White coloration %Genetic−0.81 ± 0.60−0.28 ± 0.42
      Phenotypic−0.34 ± 0.08
      SE is less than half of the correlation.
      −0.10 ± 0.06
      1 HO75 = animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree; HO50 = animals with known Holstein sire and no other known breeds in 4-generation pedigree; ME305 = mature equivalent of a 305-d lactation.
      2 SE is less than half of the correlation.
      3 NC = model did not converge.
      Table 5 presents the pedigree analysis results of the HO50 and HO75 3-trait models that included fly count. The genetic correlation between fly score and fly count was near unity in both the HO75 and HO50 data sets (0.93 and 0.95, respectively) with a high phenotypic correlation (0.74 in both data sets). From each data set, the heritability estimates for fly count (0.29 and 0.30, respectively) were higher than that of fly score (0.17 and 0.23, respectively). The high genetic correlation between fly score and fly count and limited number of records of fly count likely contributed to nonconvergence in the genomic analysis.
      Table 5Genetic parameter estimates ± SE from the three-trait pedigree models run with fly count; heritability is presented on the diagonal, genetic correlation above the diagonal and phenotypic correlation below the diagonal
      Data set
      HO75 = animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree; HO50 = animals with known Holstein sire and no other known breeds in 4-generation pedigree.
      TraitFly scoreStayabilityFly count
      HO50Fly score0.17 ± 0.040.19 ± 0.150.93 ± 0.13
      SE is less than half of the correlation.
      Stayability0.05 ± 0.03
      SE is less than half of the correlation.
      0.09 ± 0.010.42 ± 0.25
      Fly count0.74 ± 0.02
      SE is less than half of the correlation.
      0.13 ± 0.060.29 ± 0.09
      HO75Fly score0.23 ± 0.040.19 ± 0.180.95 ± 0.12
      SE is less than half of the correlation.
      Stayability0.08 ± 0.03
      SE is less than half of the correlation.
      0.07 ± 0.010.35 ± 0.28
      Fly count0.74 ± 0.03
      SE is less than half of the correlation.
      0.12 ± 0.070.30 ± 0.10
      1 HO75 = animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree; HO50 = animals with known Holstein sire and no other known breeds in 4-generation pedigree.
      2 SE is less than half of the correlation.

      Test-Day Analysis

      Table 6 summarizes the LSM of yield generated for fly scores including and excluding EBV. There was no relationship between yield and fly score when genetic merit for yield was ignored; however, a higher fly load significantly reduced test date milk and protein yield when accounting for genomic EBV. Daily milk yield declined from 32.52 kg/d to 30.24 kg/d when going from a state of very few flies (score = 0) to moderate to high levels of fly infestation (score = 3 or 4), which represented a 7% decline. Similarly, cows that received fly scores of 3 or 4 produced 0.77 kg/d of protein, which was 4.5% less contemporaries that received fly scores of 0 and produced 0.85 kg/d of protein. A numeric decline was observed in daily fat production when accounting for genomic EBV; cows that received a fly score of 0 produced 0.055 kg/d more fat than cows that scored 3 or 4. However, this difference was not statistically significant.
      Table 6Least squares means ± SE of test-day yield traits by fly score where EBV for mature equivalent of a 305-d lactation yields are included or excluded from the model
      TraitModel
      Indicates if EBV of individuals for yield traits analyzed were included in the model.
      Fly score
      Fly scores of 3 and 4 were compounded into one category notated as 3.
      P-value
      0123
      Milk (kg)EBV excluded31.86 ± 1.0731.95 ± 1.0531.86 ± 1.1031.56 ± 1.150.6353
      EBV included32.52 ± 1.0032.08 ± 0.9831.83 ± 1.0330.24 ± 1.090.0080
      P < 0.01.
      Fat (kg)EBV excluded1.06 ± 0.031.09 ± 0.021.06 ± 0.031.08 ± 0.030.9068
      EBV included1.10 ± 0.021.09 ± 0.021.06 ± 0.021.04 ± 0.030.2578
      Protein (kg)EBV excluded0.86 ± 0.030.85 ± 0.030.86 ± 0.030.84 ± 0.030.7281
      EBV included0.85 ± 0.020.83 ± 0.020.82 ± 0.020.77 ± 0.020.0011
      P < 0.01.
      1 Indicates if EBV of individuals for yield traits analyzed were included in the model.
      2 Fly scores of 3 and 4 were compounded into one category notated as 3.
      ** P < 0.01.

      GWAS

      Manhattan plots of the SNP P-values and variance over 1-Mb windows for fly score from the HO75 data are displayed in Figure 3. Table 7 summarizes the SNP identified as significant following Bonferroni correction in GWAS for fly score from each model. Each model resulted in a noticeable peak for fly score on BTA 6; genome-wide significance was detected for 4 SNP on BTA 6. The most significant SNP was located at bp 70359810 (P < 0.001). The window (69462146 bp to 70460998 bp) that explained the most fly score variance (1.7%) encompasses 3 of the significant SNP identified. The genes encompassed in this window include cysteine rich hydrophobic domain 2 (CHIC2), GS homeobox 2 (GSX2), platelet-derived growth factor receptor α (PDGFRA), and KIT, which are further described in Table 8. Additionally, 4 pseudogenes and 4 noncoding RNA (ncRNA) sequences were identified in this genomic window (Table 9).
      Figure thumbnail gr3
      Figure 3Manhattan plots generated from GWAS of fly score. (A) Fly load variance explained by 1-Mb windows in animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree (HO75). The points represent 1-Mb windows plotted by the amount of variation in fly score they explain. (B) Fly load SNP significance identified in HO75 data. The points plotted are the P-values of individual SNP effects from the HO75 data. The red line represents genome-wide significance.
      Table 7The SNP that most significantly explain fly score; P-values are presented as the equivalent significance for a single hypothesis test
      SNP IDBTALocation (bp)P-value
      Hapmap27516-BTC-0424656702466890.0433
      BovineHD06000200126703598100.00000719
      BovineHD06000200416704609980.00341
      BovineHD06000205776722595350.0139
      1 The P-values presented are Bonferroni corrected for genome-wide significance.
      Table 8Genes located in the genomic region on BTA 6 (69462146 bp to 70460998 bp) that explains the greatest proportion of variation in fly score; the same region explains the greatest proportion of variation in black and white coat coloration
      Gene symbolGene full nameLocationProtein functionReference
      CHIC2Cysteine rich hydrophobic domain 269,536,633 to 69,596,055Plasma membrane vesicle composition; associated with leukemia in humans following chromosomal translocation
      • Cools J.
      • Mentens N.
      • Marynen P.
      A new family of small, palmitoylated, membrane-associated proteins, characterized by the presence of a cysteine-rich hydrophobic motif.
      GSX2GS homeobox 269,629,306 to 69,630,906DNA-binding transcription factor that affects neural stem cell development
      • López-Juárez A.
      • Howard J.
      • Ullom K.
      • Howard L.
      • Grande A.
      • Pardo A.
      • Waclaw R.
      • Sun Y.Y.
      • Yang D.
      • Kuan C.Y.
      • Campbell K.
      • Nakafuku M.
      Gsx2 controls region-specific activation of neural stem cells and injury-induced neurogenesis in the adult subventricular zone.
      PDGFRAPlatelet-derived growth factor receptor α69,723,778 to 69,771,549Tyrosine kinase cell surface receptor for various growth factors; involvement in multiple pathways well documented
      • Lei H.
      • Kazlauskas A.
      Growth factors outside of the platelet-derived growth factor (PDGF) family employ reactive oxygen species/Src family kinases to activate PDGF receptor α and thereby promote proliferation and survival of cells.
      KIT
      Contains significant SNP identified in GWAS.
      KIT proto-oncogene, receptor tyrosine kinase70,166,681 to 70,254,049Tyrosine kinase cell surface receptor for growth factors, ATP, metal ions, and cytokines; required for melanogenesis; involvement in multiple pathways well documented
      • Miettinen M.
      • Lasota J.
      KIT (CD117): A review on expression in normal and neoplastic tissues, and mutations and their clinicopathologic correlation.
      1 Contains significant SNP identified in GWAS.
      Table 9Noncoding regions in the most significant genomic region of BTA 6 determined by GWAS
      IDLocationType
      LOC11233707969,541,388 to 69,547,064Pseudogene
      LOC11244708069,595,012 to 69,596,402ncRNA
      ncRNA = noncoding RNA.
      LOC10029697469,811,870 to 69,812,305Pseudogene
      LOC11244721169,877,227 to 69,877,333ncRNA
      LOC10029650569,912,526 to 69,912,920Pseudogene
      LOC11244720869,954,510 to 69,954,616ncRNA
      LOC11244708170,350,350 to 70,369,372ncRNA
      LOC100113856370,458,765 to 70,529,979Pseudogene
      1 ncRNA = noncoding RNA.
      The genomic regions that most explain the variation in fly score also most explain the variation in coat coloration. P-values for individual SNP and proportion of variance explained 1-Mb windows in the genome for proportion of white coat coloration are plotted in Figure 4. In total, 20 SNP were significantly associated with white coloration, 18 located on BTA 6, 1 on BTA 8, and 1 on BTA 3. The same window on BTA 6 (69462146 bp to 70460998 bp) that explained the most variation in fly score explained the most variation in proportion of white coat coloration (2.1%). All 4 SNP that were significantly associated with fly score were significantly associated with coat coloration. Their P-values are presented in Table 10.
      Figure thumbnail gr4
      Figure 4Manhattan plots generated from GWAS of proportion of white coloration. (A) White coloration variance explained by 1-Mb windows in animals with known Holstein sire and known Holstein maternal grandsire with no other known breeds in 4-generation pedigree (HO75). The points represent 1-Mb windows plotted by the amount of variation of proportion of white coloration they explain. (B) White coloration SNP significance identified in HO75 data. The points plotted are the P-values of individual SNP effects from the HO75 data. The red line represents genome-wide significance.
      Table 10The significance of SNP that most significantly explained fly score for the trait percentage of white coloration; P-values are presented as the equivalent significance for a single hypothesis test
      SNP identificationBTALocation (bp)P-value
      BovineHD06000200126703598100.0000000257
      BovineHD06000200416704609980.00000439
      BovineHD06000205776722595350.000277
      Hapmap27516-BTC-0424656702466890.000194
      1The P-values presented are Bonferroni corrected for genome-wide significance.

      Color Analyses

      Variants of Hapmap27516-BTC-042465 were considered in the following analyses because it was the only SNP significantly associated with fly score that (1) fell within the 1-Mb window that explained the most variation in scores and (2) was located on a known gene rather than a pseudogene or ncRNA sequence. The proportion of white coloration was significantly (P < 0.0001) influenced by SNP variants (Figure 5). The LSM estimates are displayed in Table 11. The proportions of coat coloration for each genotype differed significantly (P < 0.05). Animals with 0 copies of the SNP variant had the most white coloration (47%), whereas those with 2 copies of the variant had the least white coloration (12%). Cows with one copy of the SNP variant had intermediate proportions of white (15%).
      Figure thumbnail gr5
      Figure 5Boxplot of the percentage of white coat coloration in individual cattle separated by variant of Hapmap27516-BTC-042465. Pink triangles represent the mean value. Lowercase letters (a–c) represent significantly (P < 0.05) different amounts of coloration. The bottom of each box represents the 25th percentile, the thick middle line represents the median, and the top represents the 75th percentile. The whiskers represent the minimum and maximum values that fall within 150% of the interquartile range and white circles represent points that exceed such a value.
      Table 11Least squares means estimates ± SE of percentage of white coat coloration for each variant of Hapmap27516-BTC-042465
      SNP variantLSM estimateP-value
      046.95 ± 2.92<0.0001
      125.18 ± 2.82<0.0001
      211.86 ± 4.530.0095
      Table 12 displays the LS means estimates of fly score for variants of Hapmap27516-BTC-042465. Variants of Hapmap27516-BTC-042465 significantly (P < 0.05) influenced fly score, those with 0 copies of the SNP variant had the lowest fly score, whereas those with 2 copies had the highest. When proportion of white coloration was added to the model, significance for the SNP variant effect was nullified.
      Table 12Least squares means estimates ± SE of fly score by Hapmap27516-BTC-042465 variant where proportion of coloration is included or excluded from the model
      TraitModelSNP variant
      012
      Fly scoreNo color effect1.34 ± 0.19
      Values within row with different superscript are different at P < 0.05.
      Values within row with different superscript are different at P < 0.10.
      1.49 ± 0.19
      Values within row with different superscript are different at P < 0.05.
      Values within row with different superscript are different at P < 0.10.
      1.61 ± 0.21
      Values within row with different superscript are different at P < 0.05.
      White %1.43 ± 0.221.46 ± 0.211.49 ± 0.24
      a,b Values within row with different superscript are different at P < 0.05.
      c,d Values within row with different superscript are different at P < 0.10.

      DISCUSSION

      The significant rank correlation observed between fly score and fly count provided confidence that fly scoring is a valid method of quantifying horn fly loads on cows. The strong genetic correlations between the 2 traits determined from the pedigree analysis of the HO75 and HO50 further supported our decision to use fly score as the method to quantify horn fly load in our analyses. The higher heritability for fly count may suggest that it is more precise than subjective scores with relatively few categories; however, the time required to assign scores is substantially less and more amenable to assembling databases of sufficient size to enable genomic evaluations.
      Our analysis estimated a moderate heritability of fly score (0.19 to 0.25) in organic Holstein cows. Admittedly, heritability values may have been slightly underestimated because of modeling choices. Because multitrait models were used for our analysis and tertiary traits of interest were primarily continuous, linear models were used to analyze all traits though fly score was not normally distributed. Our decision to simplify modeling by using a linear model for all traits rather than modeling stayability and fly score with threshold models was justified by previous work that estimated genetic parameters of binomial hoof lesion data with both linear and threshold models (
      • Malchiodi F.
      • Koeck A.
      • Mason S.
      • Christen A.M.
      • Kelton D.F.
      • Schenkel F.S.
      • Miglior F.
      Genetic parameters for hoof health traits estimated with linear and threshold models using alternative cohorts.
      ). The EVB derived from the threshold and linear models were highly correlated (R >0.95), and heritability estimates differed by 0.01 at most (
      • Malchiodi F.
      • Koeck A.
      • Mason S.
      • Christen A.M.
      • Kelton D.F.
      • Schenkel F.S.
      • Miglior F.
      Genetic parameters for hoof health traits estimated with linear and threshold models using alternative cohorts.
      ). The results of
      • Malchiodi F.
      • Koeck A.
      • Mason S.
      • Christen A.M.
      • Kelton D.F.
      • Schenkel F.S.
      • Miglior F.
      Genetic parameters for hoof health traits estimated with linear and threshold models using alternative cohorts.
      suggest that our genetic parameter estimates would have been minimally changed if threshold modeling had been used.
      Our heritability estimate of horn fly score was greater than the heritability estimate of 0.08 in a population of Caracu beef cattle (
      • Fraga A.B.
      • Mello De Alencar M.
      • Andrade De Figueiredo L.
      • Razook A.G.
      • Noely J.
      • Cyrillo S.G.
      Genetic analysis of the infestation of females of the Caracu cattle breed by Horn Fly (Haematobia irritans irritans) (L.) (Diptera, Muscidae).
      ). An older study estimated the heritability of horn fly resistance to be between 0.59 and 0.78 based on data from a small sample of beef cattle (
      • Brown Jr., A.H.
      • Steelman C.D.
      • Johnson Z.B.
      • Rosenkrans Jr., C.F.
      • Brasuell T.M.
      Estimates of repeatability and heritability of horn fly resistance in beef cattle.
      ). The values estimated by
      • Brown Jr., A.H.
      • Steelman C.D.
      • Johnson Z.B.
      • Rosenkrans Jr., C.F.
      • Brasuell T.M.
      Estimates of repeatability and heritability of horn fly resistance in beef cattle.
      likely differed from those presented here or estimated by
      • Fraga A.B.
      • Mello De Alencar M.
      • Andrade De Figueiredo L.
      • Razook A.G.
      • Noely J.
      • Cyrillo S.G.
      Genetic analysis of the infestation of females of the Caracu cattle breed by Horn Fly (Haematobia irritans irritans) (L.) (Diptera, Muscidae).
      because of the few observations (n = 427) resulted in imprecise estimates. The heritability generated in our genomic model is substantially higher than what was estimated by
      • Fraga A.B.
      • Mello De Alencar M.
      • Andrade De Figueiredo L.
      • Razook A.G.
      • Noely J.
      • Cyrillo S.G.
      Genetic analysis of the infestation of females of the Caracu cattle breed by Horn Fly (Haematobia irritans irritans) (L.) (Diptera, Muscidae).
      , suggesting that the inclusion of genomic information allows for a more robust estimation of genetic variance over that estimated using only pedigree data.
      The small, positive, phenotypic correlation between stayability and fly score indicated that longer herd life was associated with higher fly scores. The genetic correlation between stayability and fly score was moderate but because of the high standard error we cannot infer from this data if heifers culled before maturity had an innate resistance or susceptibility to horn flies. The positive genetic correlations between fly score and ME305 milk, fat, and protein yield suggest that genetic selection for yield may have an unfavorable effect on the number of flies on a cow, though the reasons for this relationship are not clear. Volatile semiochemicals (VSC), which are pheromone-like odors that have attractant and repellant effects upon insects, can be emitted from the breath and urine of cattle and have been found to vary among individual animals (
      • Birkett M.A.
      • Agelopoulos N.
      • Jensen K.-M.V.
      • Jespersen J.B.
      • Pickett J.A.
      • Prijs H.J.
      • Thomas G.
      • Trapman J.J.
      • Wadhams L.J.
      • Woodcock C.M.
      The role of volatile semiochemicals in mediating host location and selection by nuisance and disease-transmitting cattle flies.
      ;
      • Oyarzún M.P.
      • Quiroz A.
      • Birkett M.A.
      Insecticide resistance in the horn fly: Alternative control strategies.
      ). Further, Holstein heifers emitting varying VSC also varied in horn fly load (
      • Jensen K.-M.V.
      • Jespersen J.B.
      • Birkett M.A.
      • Pickett J.A.
      • Thomas G.
      • Wadhams L.J.
      • Woodcock C.M.
      Variation in the load of the horn fly, Haematobia irritans, in cattle herds is determined by the presence or absence of individual heifers.
      ). However, the VSC identified are primarily produced by plants, fungi, and bacteria (
      • Bennett J.W.
      • Hung R.
      • Lee S.
      • Padhi S.
      ). Thus, the VSC identified as having horn fly repellant or attractant properties are likely produced by the microbiome of a cow during digestion. Digestion rates are highly correlated with milk yield (
      • Nocek J.E.
      • Russell J.B.
      Protein and energy as an integrated system. Relationship of ruminal protein and carbohydrate availability to microbial synthesis and milk production.
      ); the rate and volume of VSC emitted by high yielding cattle may be more attractive to horn flies than that of their low yield herd-mates, but this is speculative.
      The positive genetic relationship between fly load and yield may have masked the phenotypic effect of high fly load on milk yield when considering ME305 production. Our mixed model analysis of test-day records supports this theory. The inclusion of EBV in our model increased the significance of fly score for every yield trait. A numerical decrease in yield for all production traits was observed as fly score increased. The results are consistent with numerous studies that have found that milk yield decreases as horn fly infestation increases (
      • Morgan D. W.T.
      • Bailie H.D.
      A field trial to determine the effect of fly control using permethrin on milk yields in dairy cattle in the UK.
      ;
      • Jonsson N.N.
      • Mayer D.G.
      Estimation of the effects of buffalo fly (Haematobia irritans exigua) on the milk production of dairy cattle based on a meta-analysis of literature data.
      ;
      • Mays A.R.
      • Brown M.A.
      • von Tunglen D.L.
      • Rosenkrans Jr., C.F.
      Milk production traits of beef cows as affected by horn fly count and sire breed type.
      ). These data suggest that cows that are genetically inclined to higher production are more susceptible to fly infestation. However, at the time of horn fly infestation all cows in the herd experience production losses relative to their production potential. We should note that we merged phenotypes with test dates that occurred within 2 weeks of fly scoring. The effects may be even more apparent if fly scores and yield were measured on the same day.
      The genes identified within the window that explained the most genetic variation by GWAS were considered for candidate genes for fly resistance (Table 8). It should be mentioned that the 1-Mb windows of the genome plotted in Figure 3, Figure 4 were not independent from each other. The GWAS relies on linkage disequilibrium between SNP and QTL and the subsequent windows used in our analysis were also affected by linkage disequilibrium between neighboring windows. Because of this, the exact amount of variation in horn fly infestation explained by genomic regions on BTA 6 was not quantified in our GWAS. Instead, the greatest amount of linkage disequilibrium was detected in 1-Mb windows on BTA 6. Because most of the significant SNP were located in the window that explained the most variation in our analysis, we focused on genes located within that window.
      The gene CHIC2 is expressed in many bovine tissues where its protein product functions as a part of cellular vesicles and the plasma membrane (
      • Cools J.
      • Mentens N.
      • Marynen P.
      A new family of small, palmitoylated, membrane-associated proteins, characterized by the presence of a cysteine-rich hydrophobic motif.
      ;
      • Papatheodorou I.
      • Moreno P.
      • Manning J.
      • Muñoz-Pomer Fuentes A.
      • George N.
      • Fexova S.
      • Fonseca N.A.
      • Füllgrabe A.
      • Green M.
      • Huang N.
      • Huerta L.
      • Iqbal H.
      • Jianu M.
      • Mohammed S.
      • Zhao L.
      • Jarnuczak A.F.
      • Jupp S.
      • Marioni J.
      • Meyer K.
      • Petryszak R.
      • Medina C. A.P.
      • Talavera-López C.
      • Teichmann S.
      • Vizcaino J.A.
      • Brazma A.
      Expression Atlas update: From tissues to single cells.
      ). Although a simple biological relationship between this function and fly resistance cannot be made based on the current literature, future findings may implicate a functional relationship supported by this GWAS result. Similarly, the transcription factor produced by GSX2 is primarily known to function in neural development and regeneration (
      • López-Juárez A.
      • Howard J.
      • Ullom K.
      • Howard L.
      • Grande A.
      • Pardo A.
      • Waclaw R.
      • Sun Y.Y.
      • Yang D.
      • Kuan C.Y.
      • Campbell K.
      • Nakafuku M.
      Gsx2 controls region-specific activation of neural stem cells and injury-induced neurogenesis in the adult subventricular zone.
      ). Horn fly avoidance behavior, such as kicking and tail flicking, could be affected by a variant of GSX2 but there is no available literature on the variance of horn fly avoidance behavior among individuals, so this explanation is highly speculative.
      The latter genes in Table 8, PDGFRA and KIT, both encode tyrosine kinase cell surface receptors and are physically linked in multiple species (
      • Lord E.A.
      • Lumsden J.M.
      • Dodds K.G.
      • Henry H.M.
      • Crawford A.M.
      • Ansari H.A.
      • Pearce P.D.
      • Maher D.W.
      • Stone R.T.
      • Kappes S.M.
      • Beattie C.W.
      • Montgomery G.W.
      The linkage map of sheep Chromosome 6 compared with orthologous regions in other species.
      ). Both genes also have roles in the super-pathways for RET signaling, synovial fibroblast apoptosis, G-protein-coupled receptors, growth factor receptor bound protein 2 associated binder 1 signalosome, and extracellular-signal-related kinase signaling (
      • Stelzer G.
      • Rosen N.
      • Plaschkes I.
      • Zimmerman S.
      • Twik M.
      • Fishilevich S.
      • Iny Stein T.
      • Nudel R.
      • Lieder I.
      • Mazor Y.
      • Kaplan S.
      • Dahary D.
      • Warshawsky D.
      • Guan-Golan Y.
      • Kohn A.
      • Rappaport N.
      • Safran M.
      • Lancet D.
      The GeneCards Suite: From gene data mining to disease genome sequence analyses.
      ). Notably, many of these pathways control important metabolic functions (
      • Bard-Chapeau E.A.
      • Hevener A.L.
      • Long S.
      • Zhang E.E.
      • Olefsky J.M.
      • Feng G.S.
      Deletion of Gab1 in the liver leads to enhanced glucose tolerance and improved hepatic insulin action.
      ;
      • Husted A.S.
      • Trauelsen M.
      • Rudenko O.
      • Hjorth S.A.
      • Schwartz T.W.
      GPCR-mediated signaling of metabolites.
      ). Despite their roles in metabolism, it is unknown if variants of PDGFRA or KIT affect the composition or quantity of VSC emittance during metabolic processes.
      The gene KIT also plays a crucial role in melanogenesis (
      • Miettinen M.
      • Lasota J.
      KIT (CD117): A review on expression in normal and neoplastic tissues, and mutations and their clinicopathologic correlation.
      ;
      • D'Mello S. A.N.
      • Finlay G.J.
      • Baguley B.C.
      • Askarian-Amiri M.E.
      Signaling pathways in melanogenesis.
      ). It has long been considered a candidate gene for spotted coat pattern in cattle (
      • Grosz M.D.
      • MacNeil M.D.
      The “spotted” locus maps to bovine chromosome 6 in a Hereford-cross population.
      ;
      • Reinsch N.
      • Thomsen H.
      • Xu N.
      • Brink M.
      • Looft C.
      • Kalm E.
      • Brockmann G.A.
      • Grupe S.
      • Kühn C.
      • Schwerin M.
      • Leyhe B.
      • Hiendleder S.
      • Erhardt G.
      • Medjugorac I.
      • Russ I.
      • Fö M.
      • Reents R.
      • Averdunk G.
      A QTL for the degree of spotting in cattle shows synteny with the KIT locus on chromosome 6.
      ). It has since been reported to explain about 9% of the variance of the proportion of black on the coats of Holstein bulls (
      • Hayes B.J.
      • Pryce J.
      • Chamberlain A.J.
      • Bowman P.J.
      • Goddard M.E.
      Genetic architecture of complex traits and accuracy of genomic prediction: Coat colour, milk-fat percentage, and type in Holstein cattle as contrasting model traits.
      ). In some cases, horn flies have been found to prefer dark colored hosts over white hosts (
      • Franks R.E.
      • Burns E.C.
      • England N.C.
      Color preference of the horn fly, Haematobia irritans, on beef cattle.
      ;
      • Doube B.M.
      The effect of breed and coat colour on numbers of the buffalo fly Haematobia Iriitans Exigua de meijere, (diptera, muscidae) on bovine hosts.
      ). Likewise, the phenotypic correlations we observed between fly score and coat color of organic Holstein cows found that cows with a greater proportion of white in their coat had lower horn fly scores than their herd-mates that had greater proportions of black coloration. Red-coated cattle were observed to host fewer flies than black cattle in one study and fewer flies than bay and orange cattle in another (
      • Schreiber E.T.
      • Campbell J.B.
      Horn fly (Diptera: Muscidae) distribution on cattle as influenced by host color and time of day.
      ;
      • Fraga A.B.
      • Mello De Alencar M.
      • Andrade De Figueiredo L.
      • Razook A.G.
      • Noely J.
      • Cyrillo S.G.
      Genetic analysis of the infestation of females of the Caracu cattle breed by Horn Fly (Haematobia irritans irritans) (L.) (Diptera, Muscidae).
      ). Contrarily, no significance differences were observed in fly load between red and black beef cattle or between solid black Jersey-Holstein crossbreds and spotted Holstein steers (
      • Franks R.E.
      • Burns E.C.
      • England N.C.
      Color preference of the horn fly, Haematobia irritans, on beef cattle.
      ;
      • Guglielmone A.A.
      • Volpogni M.M.
      • Castro H.
      • Mangold A.J.
      • Anziani O.S.
      A study of relative horn fly, Haematobia irritans (Diptera: Muscidae), abundance on Holstein steers and steers of two Holstein crosses.
      ). The latter may be due to breed and coat color serving as confounding variables. Most recently, black cattle painted with black and white stripes attracted fewer biting flies than their counterparts painted only with black stripes or with no paint at all, implying that color or patterning may affect biting fly attractiveness (
      • Kojima T.
      • Oishi K.
      • Matsubara Y.
      • Uchiyama Y.
      • Fukushima Y.
      • Aoki N.
      • Sato S.
      • Masuda T.
      • Ueda J.
      • Hirooka H.
      • Kino K.
      Cows painted with zebra-like striping can avoid biting fly attack.
      ).
      Our analysis of a significant SNP identified, Hapmap27516-BTC-042465, which falls on KIT, further supports its candidacy as a gene that plays a factor in coat coloration and fly resistance in Holstein cattle. Cows with 0 copies of the SNP variant had fewer horn flies and more white coat coloration, whereas those with 2 copies had more flies and more black coat coloration. The loss of SNP variant significance in models of fly score that included coat coloration effects suggest that the variation explained by SNP variant was explained similarly by coat coloration.
      The results presented here show that there is a relationship between horn fly infestation in cattle and coat color. It appears that Holsteins with more white coloration are infested with fewer horn flies whereas those with more black coloration are infested with more. This is supported by the strong negative phenotypic correlation between white coloration and fly score. Further, KIT, which is known to be associated with coloration, is implicated as a genetic driver of both coat color and horn fly infestation. More coat color phenotypes need to be collected from red and white Holsteins and other breeds of cattle to determine if KIT influences the coloration and subsequent horn fly loads of cattle with different colored hides. Both coat coloration and horn fly resistance are likely polygenic traits; thus, more phenotypes of fly load could be beneficial to identify other candidate genes that may influence horn fly resistance by mechanisms in addition to coloration.
      Finally, it is worth mentioning that the pseudogenes and ncRNA sequences identified in Table 9 could also explain variation in fly load. Although both have long been considered “junk” components of the DNA sequence, recent works have revealed that they play important roles in gene regulation (
      • Tutar Y.
      Psuedogenes.
      ;
      • Palazzo A.F.
      • Lee E.S.
      Non-coding RNA: What is functional and what is junk?.
      ). Investigation into these effects are beyond the scope of this study but their potential to influence genetic variation is notable.
      Horn fly resistance is a moderately heritable trait in Holstein cattle that could be influenced by selection. Fly scoring is a fast and accurate method of quantifying horn fly load that could be used to develop a database of sufficient size for genomic analysis. There was a positive and unfavorable genetic correlation between fly load and ME305 milk and fat yields, which indicates continued selection for higher yield may increase the fly susceptibility of pastured cattle over time. When accounting for breeding values, test-day yields of milk and protein were lower in cattle highly infested with flies. Fly infestation appeared to be a largely polygenic trait but the genomic region between 69462146 bp and 70460998 bp was significantly associated with horn fly resistance in Holstein cattle. Coat coloration appears to influence susceptibility to horn fly infestation; Holstein cows with greater proportions of white coloration attracted fewer horn flies than those that were mostly black. The relationship between horn fly infestation and coat coloration supports KIT as a candidate gene for horn fly resistance.
      In the future, sequencing of cattle consistently highly and lowly infested with horn flies could be used to identify any causative mutations that may be located on KIT or other genes. The incorporation of fly scores into data management records could be implemented so that fly resistance could be developed into a selection index trait. Horn fly resistant cattle would especially benefit organic and pasture-based operations by preventing production losses associated with infestation and reducing disease spread by flies. Both organic and conventional producers could reduce costs and labor associated with horn fly control by selecting for horn fly resistance in their herds. Fly resistance further has positive implications toward animal welfare and products from fly resistant cattle could be more desirable to the consumer.

      ACKNOWLEDGMENTS

      The authors thank the organic dairy herds that participated in this study. Charlotte Natali, Sarah Kunkel, George Demers, and Abigail Jantzi (Penn State University, University Park, PA) are recognized for their work in counting flies and scoring coat colors. This work was supported by the USDA National Institute of Food and Agriculture Organic Agriculture Research and Extension Initiative competitive grant number 2016-51300-25862. The authors have not stated any conflicts of interest.

      REFERENCES

        • Aguilar I.
        • Tsuruta S.
        • Legarra A.
        PREGSF90-POSTGSF90: Computational tools for the implementation of single-step genomic selection and genome-wide association with ungenotyped individuals in BLUPF90 programs. Proc. 10th World Congr. Genet. Appl. Livest. Prod. American Society of Animal Science.
        • Árnason T.
        • Albertsdóttir E.
        • Fikse W.F.
        • Eriksson S.
        • Sigurdsson Á.
        Estimation of genetic parameters and response to selection for a continuous trait subject to culling before testing.
        J. Anim. Breed. Genet. 2012; 129 (22225584): 50-59
        • Bard-Chapeau E.A.
        • Hevener A.L.
        • Long S.
        • Zhang E.E.
        • Olefsky J.M.
        • Feng G.S.
        Deletion of Gab1 in the liver leads to enhanced glucose tolerance and improved hepatic insulin action.
        Nat. Med. 2005; 11 (15821749): 567-571
        • Bennett J.W.
        • Hung R.
        • Lee S.
        • Padhi S.
        Hock B. Fungal and Bacterial Volatile Organic Compounds: An Overview and Their Role. 2nd ed. Springer, 2012
        • Birkett M.A.
        • Agelopoulos N.
        • Jensen K.-M.V.
        • Jespersen J.B.
        • Pickett J.A.
        • Prijs H.J.
        • Thomas G.
        • Trapman J.J.
        • Wadhams L.J.
        • Woodcock C.M.
        The role of volatile semiochemicals in mediating host location and selection by nuisance and disease-transmitting cattle flies.
        Med. Vet. Entomol. 2004; 18 (15641996): 313-322
        • Brock C.C.
        • Pempek J.A.
        • Jackson-Smith D.
        • Weaver K.
        • da Costa L.
        • Habing G.G.
        Organic dairy producer experiences and decisions related to disease prevention and treatment.
        J. Dairy Sci. 2021; 104 (33516550): 5867-5880
        • Brown Jr., A.H.
        • Steelman C.D.
        • Johnson Z.B.
        • Rosenkrans Jr., C.F.
        • Brasuell T.M.
        Estimates of repeatability and heritability of horn fly resistance in beef cattle.
        J. Anim. Sci. 1992; 70 (1526906): 1375-1381
        • Byford R.L.
        • Craig M.E.
        • Crosby B.L.
        A review of ectoparasites and their effect on cattle production.
        J. Anim. Sci. 1992; 70 (1347767): 597-602
        • Clark K.
        • Karsch-Mizrachi I.
        • Lipman D.J.
        • Ostell J.
        • Sayers E.W.
        GenBank.
        Nucleic Acids Res. 2016; 44 (26590407): D67-D72
        • Cole J.B.
        • VanRaden P.M.
        A Manual for Use of BESTPRED: A Program for Estimation of Lactation Yield and Persistency Using Best Prediction.
        USDA AIPL, 2009
        • Cools J.
        • Mentens N.
        • Marynen P.
        A new family of small, palmitoylated, membrane-associated proteins, characterized by the presence of a cysteine-rich hydrophobic motif.
        FEBS Lett. 2001; 492 (11257495): 204-209
        • D'Mello S. A.N.
        • Finlay G.J.
        • Baguley B.C.
        • Askarian-Amiri M.E.
        Signaling pathways in melanogenesis.
        Int. J. Mol. Sci. 2016; 17 (27428965)1144
        • Doube B.M.
        The effect of breed and coat colour on numbers of the buffalo fly Haematobia Iriitans Exigua de meijere, (diptera, muscidae) on bovine hosts.
        Aust. J. Entomol. 1984; 23: 39-45
        • Falconer D.S.
        • Mackay T. F.C.
        Introduction to Quantitative Genetics.
        4th ed. Prentice Hall, 1996
        • Foil L.D.
        • Hogsette J.A.
        Biology and control of tabanids, stable flies and horn flies.
        Rev. Sci. Tech. 1994; 13: 1125-1158
        • Fraga A.B.
        • Mello De Alencar M.
        • Andrade De Figueiredo L.
        • Razook A.G.
        • Noely J.
        • Cyrillo S.G.
        Genetic analysis of the infestation of females of the Caracu cattle breed by Horn Fly (Haematobia irritans irritans) (L.) (Diptera, Muscidae).
        Genet. Mol. Biol. 2005; 28: 242-247
        • Franks R.E.
        • Burns E.C.
        • England N.C.
        Color preference of the horn fly, Haematobia irritans, on beef cattle.
        J. Econ. Entomol. 1964; 57: 371-372
        • Gilmour A.R.
        • Gogel B.J.
        • Cullis B.R.
        • Welham S.J.
        • Thompson R.
        ASReml User Guide Release 4.1.
        VSN International Ltd, 2015
        • Grosz M.D.
        • MacNeil M.D.
        The “spotted” locus maps to bovine chromosome 6 in a Hereford-cross population.
        J. Hered. 1999; 90 (9987932): 233-236
        • Guglielmone A.A.
        • Volpogni M.M.
        • Castro H.
        • Mangold A.J.
        • Anziani O.S.
        A study of relative horn fly, Haematobia irritans (Diptera: Muscidae), abundance on Holstein steers and steers of two Holstein crosses.
        Vet. Parasitol. 2002; 109 (12383633): 141-145
        • Hardie L.C.
        • Heins B.J.
        • Dechow C.D.
        Genetic parameters for stayability of Holsteins in US organic herds.
        J. Dairy Sci. 2021; 104 (33589261): 4507-4515
        • Hayes B.J.
        • Pryce J.
        • Chamberlain A.J.
        • Bowman P.J.
        • Goddard M.E.
        Genetic architecture of complex traits and accuracy of genomic prediction: Coat colour, milk-fat percentage, and type in Holstein cattle as contrasting model traits.
        PLoS Genet. 2010; 6 (20927186)e1001139
        • Houle D.
        • Meyer K.
        Estimating sampling error of evolutionary statistics based on genetic covariance matrices using maximum likelihood.
        J. Evol. Biol. 2015; 28 (26079756): 1542-1549
        • Hudson G. F.S.
        • Van Vleck L.D.
        Relationship between production and stayability in Holstein cattle.
        J. Dairy Sci. 1981; 64: 2246-2250
        • Jensen K.-M.V.
        • Jespersen J.B.
        • Birkett M.A.
        • Pickett J.A.
        • Thomas G.
        • Wadhams L.J.
        • Woodcock C.M.
        Variation in the load of the horn fly, Haematobia irritans, in cattle herds is determined by the presence or absence of individual heifers.
        Med. Vet. Entomol. 2004; 18 (15347395): 275-280
        • Jonsson N.N.
        • Mayer D.G.
        Estimation of the effects of buffalo fly (Haematobia irritans exigua) on the milk production of dairy cattle based on a meta-analysis of literature data.
        Med. Vet. Entomol. 1999; 13 (10608225): 372-376
        • Kerlin R.L.
        • Allingham P.G.
        Acquired immune response of cattle exposed to buffalo fly (Haematobia irritans exigua).
        Vet. Parasitol. 1992; 43 (1496794): 115-129
        • Kojima T.
        • Oishi K.
        • Matsubara Y.
        • Uchiyama Y.
        • Fukushima Y.
        • Aoki N.
        • Sato S.
        • Masuda T.
        • Ueda J.
        • Hirooka H.
        • Kino K.
        Cows painted with zebra-like striping can avoid biting fly attack.
        PLoS One. 2019; 14 (31581218)e0223447
        • Husted A.S.
        • Trauelsen M.
        • Rudenko O.
        • Hjorth S.A.
        • Schwartz T.W.
        GPCR-mediated signaling of metabolites.
        Cell Metab. 2017; 25: 777-796
        • Lei H.
        • Kazlauskas A.
        Growth factors outside of the platelet-derived growth factor (PDGF) family employ reactive oxygen species/Src family kinases to activate PDGF receptor α and thereby promote proliferation and survival of cells.
        J. Biol. Chem. 2009; 284 (19126548): 6329-6336
        • López-Juárez A.
        • Howard J.
        • Ullom K.
        • Howard L.
        • Grande A.
        • Pardo A.
        • Waclaw R.
        • Sun Y.Y.
        • Yang D.
        • Kuan C.Y.
        • Campbell K.
        • Nakafuku M.
        Gsx2 controls region-specific activation of neural stem cells and injury-induced neurogenesis in the adult subventricular zone.
        Genes Dev. 2013; 27 (23723414): 1272-1287
        • Lord E.A.
        • Lumsden J.M.
        • Dodds K.G.
        • Henry H.M.
        • Crawford A.M.
        • Ansari H.A.
        • Pearce P.D.
        • Maher D.W.
        • Stone R.T.
        • Kappes S.M.
        • Beattie C.W.
        • Montgomery G.W.
        The linkage map of sheep Chromosome 6 compared with orthologous regions in other species.
        Mamm. Genome. 1996; 7 (8661727): 373-376
        • Lysyk T.J.
        Effect of larval rearing temperature and maternal photoperiod on diapause in the horn fly (Diptera: Muscidae).
        Environ. Entomol. 1992; 21: 1134-1138
        • Malchiodi F.
        • Koeck A.
        • Mason S.
        • Christen A.M.
        • Kelton D.F.
        • Schenkel F.S.
        • Miglior F.
        Genetic parameters for hoof health traits estimated with linear and threshold models using alternative cohorts.
        J. Dairy Sci. 2017; 100 (28131577): 2828-2836
        • Mays A.R.
        • Brown M.A.
        • von Tunglen D.L.
        • Rosenkrans Jr., C.F.
        Milk production traits of beef cows as affected by horn fly count and sire breed type.
        J. Anim. Sci. 2014; 92 (24492544): 1208-1212
        • McKay L.
        • DeLong K.L.
        • Schexnayder S.
        • Griffith A.P.
        • Taylor D.B.
        • Olafson P.
        • Trout Fryxell R.T.
        Cow-calf producers' willingness to pay for bulls resistant to horn flies (Diptera: Muscidae).
        J. Econ. Entomol. 2019; 112 (30768165): 1476-1484
        • Miettinen M.
        • Lasota J.
        KIT (CD117): A review on expression in normal and neoplastic tissues, and mutations and their clinicopathologic correlation.
        Appl. Immunohistochem. Mol. Morphol. 2005; 13 (16082245): 205-220
        • Misztal I.
        • Tsuruta S.
        • Lourenco D.
        • Aguilar I.
        • Legarra A.
        • Vitezica Z.
        Manual for BLUPF90 family of programs.
        http://nce.ads.uga.edu/wiki/doku.php?id=documentation
        Date: 2018
        Date accessed: January 29, 2021
        • Moon R.D.
        Chapter 17: Muscid flies (Muscidae).
        in: Mullen G.R. Durden L.A. Medical and Veterinary Entomology. 3rd ed. Academic Press, 2019: 345-364
        • Morgan D. W.T.
        • Bailie H.D.
        A field trial to determine the effect of fly control using permethrin on milk yields in dairy cattle in the UK.
        Vet. Rec. 1980; 106 (7445328): 121-123
        • Nocek J.E.
        • Russell J.B.
        Protein and energy as an integrated system. Relationship of ruminal protein and carbohydrate availability to microbial synthesis and milk production.
        J. Dairy Sci. 1988; 71: 2070-2107
        • Norman H.D.
        • VanRaden P.M.
        • Megonigal J.H.
        • Dürr J.W.
        • Cooper T.A.
        Breed base representation in dairy animals of five breeds.
        J. Anim. Sci. 2016; 94: 155-156
        • O'Leary N.A.
        • Wright M.W.
        • Brister J.R.
        • Ciufo S.
        • Haddad D.
        • McVeigh R.
        • Rajput B.
        • Robbertse B.
        • Smith-White B.
        • Ako-Adjei D.
        • Astashyn A.
        • Badretdin A.
        • Bao Y.
        • Blinkova O.
        • Brover V.
        • Chetvernin V.
        • Choi J.
        • Cox E.
        • Ermolaeva O.
        • Farrell C.M.
        • Goldfarb T.
        • Gupta T.
        • Haft D.
        • Hatcher E.
        • Hlavina W.
        • Joardar V.S.
        • Kodali V.K.
        • Li W.
        • Maglott D.
        • Masterson P.
        • McGarvey K.M.
        • Murphy M.R.
        • O'Neill K.
        • Pujar S.
        • Rangwala S.H.
        • Rausch D.
        • Riddick L.D.
        • Schoch C.
        • Shkeda A.
        • Storz S.S.
        • Sun H.
        • Thibaud-Nissen F.
        • Tolstoy I.
        • Tully R.E.
        • Vatsan A.R.
        • Wallin C.
        • Webb D.
        • Wu W.
        • Landrum M.J.
        • Kimchi A.
        • Tatusova T.
        • DiCuccio M.
        • Kitts P.
        • Murphy T.D.
        • Pruitt K.D.
        Reference sequence (RefSeq) database at NCBI: Current status, taxonomic expansion, and functional annotation.
        Nucleic Acids Res. 2016; 44 (26553804): D733-D745
        • Oyarzún M.P.
        • Quiroz A.
        • Birkett M.A.
        Insecticide resistance in the horn fly: Alternative control strategies.
        Med. Vet. Entomol. 2008; 22 (18816268): 188-202
        • Palazzo A.F.
        • Lee E.S.
        Non-coding RNA: What is functional and what is junk?.
        Front. Genet. 2015; 6 (25674102): 2
        • Papatheodorou I.
        • Moreno P.
        • Manning J.
        • Muñoz-Pomer Fuentes A.
        • George N.
        • Fexova S.
        • Fonseca N.A.
        • Füllgrabe A.
        • Green M.
        • Huang N.
        • Huerta L.
        • Iqbal H.
        • Jianu M.
        • Mohammed S.
        • Zhao L.
        • Jarnuczak A.F.
        • Jupp S.
        • Marioni J.
        • Meyer K.
        • Petryszak R.
        • Medina C. A.P.
        • Talavera-López C.
        • Teichmann S.
        • Vizcaino J.A.
        • Brazma A.
        Expression Atlas update: From tissues to single cells.
        Nucleic Acids Res. 2020; 48 (31665515): D77-D83
        • Reinsch N.
        • Thomsen H.
        • Xu N.
        • Brink M.
        • Looft C.
        • Kalm E.
        • Brockmann G.A.
        • Grupe S.
        • Kühn C.
        • Schwerin M.
        • Leyhe B.
        • Hiendleder S.
        • Erhardt G.
        • Medjugorac I.
        • Russ I.
        • Fö M.
        • Reents R.
        • Averdunk G.
        A QTL for the degree of spotting in cattle shows synteny with the KIT locus on chromosome 6.
        J. Hered. 1999; 90 (10589513): 629-634
        • Schneider C.A.
        • Rasband W.S.
        • Eliceiri K.W.
        NIH image to ImageJ: 25 years of image analysis.
        Nat. Methods. 2012; 9 (22930834): 671-675
        • Schreiber E.T.
        • Campbell J.B.
        Horn fly (Diptera: Muscidae) distribution on cattle as influenced by host color and time of day.
        Environ. Entomol. 1986; 15: 1307-1309
        • Stelzer G.
        • Rosen N.
        • Plaschkes I.
        • Zimmerman S.
        • Twik M.
        • Fishilevich S.
        • Iny Stein T.
        • Nudel R.
        • Lieder I.
        • Mazor Y.
        • Kaplan S.
        • Dahary D.
        • Warshawsky D.
        • Guan-Golan Y.
        • Kohn A.
        • Rappaport N.
        • Safran M.
        • Lancet D.
        The GeneCards Suite: From gene data mining to disease genome sequence analyses.
        Curr. Protoc. Bioinforma. 2016; 2016: 1.30.1-1.30.33
        • Tutar Y.
        Psuedogenes.
        Comp. Funct. Genomics. 2012; 2012 (22611337)424526
        • USDA-NOP (United States Department of Agriculture National Organic Program)
        CFR. § 205.
        • VanRaden P.M.
        • Cooper T.A.
        Genomic evaluations and breed composition for crossbred U.S. dairy cattle.
        Interbull Bull. 2015; 49: 19-23
        • VanRaden P.M.
        • O'Connell J.R.
        • Wiggans G.R.
        • Weigel K.A.
        Genomic evaluations with many more genotypes.
        Genet. Sel. Evol. 2011; 43 (21366914): 10
        • Wiggans G.R.
        • VanRaden P.M.
        • Null D.J.
        • Cole J.B.
        Genomic predictions using more markers and gene tests.
        J. Dairy Sci. 2019; 102 (Abstr.): 397

      CHORUS Manuscript

      View Open Manuscript