Identification of loci associated with susceptibility to Mycobacterium avium subspecies paratuberculosis (Map) tissue infection in cattle Open Access

Abstract

Johne's disease is a contagious bacterial infection of cattle caused by Mycobacterium avium ssp. paratuberculosis (Map). A previous genome-wide association analysis (GWAA) in Holstein cattle identified QTL on BTA3 and BTA9 that were highly associated (P < 5 · 10^-7) and on BTA1, BTA16, and BTA21 that were moderately associated (P < 5 · 10^-5) with Map tissue infection. The objectives of this study were to validate previous GWAA results in Jersey cattle (n = 57), Holstein cattle from the Pacific Northwest (PNW, n = 205) and a combined Holstein population from the PNW and the Northeast (PNW + NE, n = 423), and also identify new loci associated with Map tissue infection. DNA was genotyped using the Illumina BovineSNP50 BeadChip, and the PNW + NE data was also imputed to whole genome sequence level using Run4 of the 1000 Bull Genomes project with Beagle v 4.1 and FImpute. Cases were ileocecal node positive and controls were negative for Map by quantitative PCR (qPCR). Individuals were removed for SNP call rate < 90%, and SNP were removed for genotype call rate < 90% or minor allele frequency < 1%. For the Jersey, PNW, and PNW + NE, GWAA were conducted using an allelic dosage model. For the PNW and the PNW + NE, an additional efficient mixed-model association eXpedited (EMMAX) analysis was performed using additive, dominance and recessive models. Seven QTL on BTA22 were identified in the Jersey population with the most significant (P = 4.45 · 10^-6) located at 21.7 megabases (Mb). Six QTL were associated in the PNW and the PNW + NE analyses, including a QTL previously identified on BTA16 in the NE population. The most significant locus for the PNW was located on BTA21 at 61 Mb (P = 8.61 · 10^-8) while the most significant locus for the PNW + NE was on BTA12 at 90 Mb (P = 2.33 · 10^-5). No additional QTL were identified with the imputed GWAA. Putative positional candidate genes were identified within 50 kb 5' and 3' of each QTL. Two positional candidate genes were identified in Jersey cattle, 1 identified in the PNW and 8 in the PNW + NE populations. Many identified positional candidate genes are involved in signal transduction, have immunological functions, or have putative functional relevance in Map entry into host cells. This study supported 2 previously identified SNP within a QTL on BTA16 and identified 16 new QTL, including 2 found in the PNW and the PNW+NE, associated with Map tissue infection.

INTRODUCTION

Johne's disease, also known as bovine paratuberculosis, is an infectious and incurable disease caused by the bacterium Mycobacterium avium ssp. paratuberculosis (Map) that commonly affects ruminant species such as cattle, goats, and deer. Symptoms of the disease in cattle range from decreased milk production in subclinical cases to severe diarrhea, progressive emaciation and death in clinically infected cattle. In US dairy cattle, Map infection has been associated with decreased production resulting in estimated annual losses exceeding $1.5 billion (Harris and Barletta, 2001). No effective vaccination protocols for Johne's disease currently exist, instead voluntary management and prevention programs help producers prevent the spread of infection within and between herds (Geraghty et al., 2014). Current testing methods have limited sensitivity for detecting cattle in the early stages of infection, a major contributor to the lack of progress that control programs have had in lowering disease prevalence and infection rates within herds [World Organization for Animal Health (OIE), 2013]. An alternative approach to limiting the spread of the disease is needed. Several studies have identified loci associated with susceptibility to Map infection, though few markers have been validated between studies. The objectives of this study were to attempt to validate previously identified loci and identify new loci associated with risk of Map tissue infection in 3, independent cattle populations using the Illumina (San Diego, CA) BovineSNP50 BeadChip assay and imputed whole genome sequence (WGS) data. The populations consisted of Holstein and Jersey cattle, and individuals were classified as either Map tissue positive or Map tissue negative. A case-control study was conducted to identify loci that were associated with Map tissue infection within dairy cattle which may be used in the future to aid Johne's disease control programs through marker assisted selection.

MATERIALS AND METHODS

All animal care and sample collections were approved and performed in accordance with Washington State University Institutional Animal Care and Use Committee (study # 04073). Tissue samples were taken postmortem.

Study Populations

Three populations of cattle were used in this project: a Jersey population and 2 Holstein populations, 1 from the Pacific Northwest (PNW) and 1 from the Northeast (NE). A total of 57 Jersey cattle were used in the Jersey genome-wide association analysis (GWAA) including 9 steers from Pennsylvania and 48 cull cows from Oregon. The mean age of the control cattle was 52.6 mo and the mean age of Map positive cattle was 46.2 mo. The mean age of the cases and controls did not differ (P = 0.48) among the Jersey cattle. The PNW Holstein population consisted of 205 cull cows from several southern Idaho dairies, whose ages were unknown. The NE Holstein population consisted of 245 cows from dairies in New York, Pennsylvania, and Vermont as previously described (Settles et al., 2009). The mean age of the control cattle was 58.5 mo, while the mean age of the case cattle was 60.7 mo. The mean ages of the cases and controls did not differ (P = 0.44) among the NE Holstein cattle.

Cattle within the Jersey and the PNW Holstein populations were tested for the presence of Map using quantitative PCR (qPCR) following Map isolation from ileocecal lymph nodes. Mycobacterium avium ssp. paratuberculosis DNA was isolated from tissue samples following the protocol for the MagMax Total Nucleic Acid Isolation Kit (Ambion, Grand Island, NY) and then qPCR for Map were performed using the AgPath-ID MAP Reagent Kit (Applied Biosystems, Grand Island, NY). AgPath uses a Map-specific PCR primer set to determine if Map is present within a sample. The presence of Map is measured using a threshold cycle (C_T) value where samples that are amplified before the CT value have higher amounts of Map and samples that do not amplify indicate that no Map is present. Using this method, individuals were considered Map positive (cases) if the CT values were < 37. Individuals were classified as controls if no amplification was seen. All samples were compared to a standard curve generated by spiked samples of known quantities of Map for each qPCR. Using this method, the sensitivity of the qPCR method ranged between 90 and 95% while the specificity ranged between 99 and 100%. A randomly subsampled set of the cases (n = 5) and controls (n = 5) were also sent to an American Association of Veterinary Laboratory Diagnosticians accredited National Animal Health Laboratory for validation. There were no discordances between Map tissue positive and negative samples identified in the research laboratory and the diagnostic laboratory. After qPCR testing, 16 Jersey cattle were classified as cases and 41 as controls. For the PNW Holsteins, 75 cows were identified as cases and 130 cows were classified as controls.

Tissue infection of the NE Holstein cattle was determined by the culturing of Map from tissue as previously described (Whitlock et al., 2000). Individuals were considered cases if any cultured sample had > 0 CFU/g of tissue; individuals with 0 CFU/g in all cultured samples were considered controls. Ninety-four NE Holsteins were identified as cases and 138 were identified as controls. Thirteen cattle were not tissue cultured, were classified as untested and were not used in the study.

Genotyping

Bovine DNA was extracted from the ileum or ileocecal node using the Puregene DNA extraction kit following the manufacturer's instructions (Qiagen, Germantown, MD). The DNA was quantified using a NanoDrop 1000 spectrophotometer (ThermoScientific, Wilmington, DE). All cattle were genotyped at GeneSeek (a Neogen company, Lincoln, NE) using the Illumina BovineSNP50 v1 (NE cattle) or v2 (Jersey and PNW cattle) BeadChip (San Diego, CA) that contain approximately 54,009 SNP. The genomic coordinates of each of the SNP and the mapping of A/B alleles to reference/non-reference assembly alleles for the v2 assay were based on the forward strand of the UMD3.1 reference genome (ftp://ftp.cbcb.umd.edu/pub/data/Bos_taurus/).

Two meta-analyses using Holsteins from the PNW (n = 178) and the NE (n = 245) populations were conducted. One GWAA meta-analysis was performed using genotypes from the Illumina BovineSNP50 BeadChip and the second meta-analysis was performed from genotypes imputed to WGS level. To allow for a meta-analysis of the PNW and NE populations, the SNP manifests for both versions of the BeadChip were intersected and SNP that were not common to both assays were removed. This resulted in a total of 52,815 SNP being available for analysis in the combined PNW and NE Holstein GWAA meta-analysis using the BovineSNP50 data.

Imputation

The combined PNW and NE Holstein BovineSNP50 data were imputed to WGS variation using a 2-stage process (van Binsbergen et al., 2014; Taylor et al., 2016). The BovineSNP50 genotypes, filtered based on Holstein population call rate of > 90% and minor allele frequency > 1%, were converted from A/B allele calls provided by Illumina to a reference versus alternate allelic dosage format. This format was based on a corrected Illumina manifest for the assay based on a set of cattle that had previously been genotyped with the BovineSNP50 assay and that had also been sequenced. Using Beagle v4.1 (Browning and Browning, 2016), the newly formatted genotypes were phased and imputed to the density of the Illumina BovineHD BeadChip (San Diego, CA), which contains 777,962 SNP. This was done using a group of 2703 Holstein cattle that had previously been genotyped on the BovineHD BeadChip as a reference (Neibergs et al., 2014). The imputed BovineHD data were then imputed to WGS level (35,431,201 indels and SNP) with FImpute using phased Run4 WGS data for 1,147 previously sequenced cattle from the 1,000 Bull Genomes Project as a reference (Sargolzaei et al., 2014).

Quality Control

Jersey Population

Prior to analysis, genotype data were quality control filtered. Cattle were first filtered for genotype call rate and 4 cattle with > 10% of their genotypes missing were removed. Two cattle were genotyped twice to assess genotyping accuracy which was 99.9%. One pair of samples with duplicated genotypes had both samples removed due to discrepancy between the reported phenotypes between members of the pair. The second pair of samples with duplicated genotypes had concordant phenotypes, so only a single duplicate was removed from the analysis as these cattle were inadvertently genotyped twice.

A principal component analysis (PCA) was used to determine if population stratification existed within the Jersey population. The PCA identified 3 groupings of cattle that were separated from the main cluster of cattle (data not shown). These cattle were identified as half-sibs when an identity by descent (IBD) matrix was computed. A genomic relationship matrix (GRM) was generated using the efficient mixed-model association eXpedited (EMMAX) software and the genomic inflation factor (λGC) was estimated to determine if the model was able to correct for the pedigree-based stratification identified by the PCA. No cattle required removal from the EMMAX-GRM analyses because λ_GC = 1.0 indicating that the significance values for the test statistics were not inflated by population stratification. A genotype-based sex check analysis was also performed to confirm or refute the phenotypic designation of the 9 males and 48 females. Cattle were identified as genetic females if the heterozygosity rate of the X chromosome was greater than 0.2 and were identified as genetic males if the X chromosome heterozygosity rate was less than 0.2 (the X chromosome is significantly misassembled in UMD3.1 and is known to include autosomal loci; Turner et al., 2011). The genetic and anatomical designations of sex were concordant in all cattle, so none were removed. After quality control filtering, a total of 50 cattle (15 cases and 35 controls) remained for the Jersey GWAA. A chi-square P-value was calculated to determine if location (Pennsylvania or Oregon) was a confounding factor that needed to be included as a covariate in the mixed models employed for the GWAA. There was sufficient evidence (P = 0.01) to support the hypothesis that location contributed to risk of Map infection in the Jersey cattle. Quality control was also performed for SNP and 1,434 SNP with a call rate < 90% were removed along with 11,736 SNP with minor allele frequencies (MAF) < 1%. All SNP passed the Hardy-Weinberg equilibrium (HWE) test with a significance threshold of P < 1 · 10^-50.

Pacific Northwest Holstein Population

Quality control was performed on both cattle and SNP prior to the GWAA in the PNW Holstein population in a manner that was similar to that described for the Jersey population. Nine cattle were removed for having a genotyping call rate of < 90%. Three cattle were genotyped twice and their replicated genotypes were 99.9% concordant. Of the 3 pairs of samples with duplicate genotypes, 2 pairs had both samples removed due to the reported phenotypes being discrepant between both members of each pair, while the third pair had concordant phenotypes and so a single duplicate was removed. The duplicated samples were due to sending the same samples twice for genotyping. Principal component analysis identified 13 cattle that distinctly grouped away from the majority of the other cattle. A GRM was created to examine the extent of relatedness among these outliers, but none of the cattle were highly related. As for the Jersey population, λ_GC was estimated and determined to be 1.0, so no further removal of cattle to remove population stratification effects was warranted. A sex check analysis revealed that 3 of the Holstein cattle were identified as having different anatomical and genetic sexual designations and were removed from the data. After filtering, 188 cattle (70 cases and 118 controls) remained for the GWAA. Quality control filtering removed 646 SNP with a genotype call rate < 90% and 6,877 SNP with a MAF < 1%. All SNP passed the HWE test threshold of P < 1 · 10^-50.

Combined Pacific Northwest and Northeast Holstein Meta-Analyses

Twenty-seven cattle were removed prior to the meta-analyses for failing cattle quality control measures in the individual (PNW and NE) GWAA analyses. The initial meta-analysis was based on the Illumina BovineSNP50 BeadChip genotypes. As cattle that had previously failed to meet quality control metrics in the individual Holstein GWAA were removed, no further cattle were removed for quality control (n = 409; 162 cases and 247 controls). The PCA identified 12 cattle that were distinctly grouped away from the majority of the cattle. These cattle were identified as half-sibs using a GRM, but were not removed from the analysis as the λ_GC was calculated to be 1.0 indicating that significance values were not overestimated due to issues of population stratification. After quality control, 409 cattle (162 cases and 247 controls) remained for analysis. A chi-square P-value was calculated to determine if location (NE or PNW) influenced risk of Map infection, however, there was no evidence (P = 0.277) of a regional effect. After evaluation of the genotype call rates, 503 SNP that had a genotype call rate < 90% and 6,953 SNP with a MAF < 1% were removed from the GWAA. Three SNP were removed for failing the HWE test threshold of P < 1 · 10^-50.

Genotypes from the Illumina BovineSNP50 BeadChip were subsequently imputed to WGS level for the cattle that remained after the initial quality control filtering (n = 409; 162 cases and 247 controls) and were used in the second meta-analysis. Quality control for SNP in the combined population resulted in removal of 3,995 SNP with a genotype call rate < 90% and 19,650,902 SNP with a MAF < 1%. A HWE test was performed using a HWE significance threshold of P < 1 · 10^-50 and a total of 24,253 SNP were removed from further analysis.

Statistical Analysis

All analyses for the Jersey, PNW Holstein and the combined PNW and NE Holstein populations were performed using EMMAX-GRM executed in the SNP & Variation Suite v8 environment (Golden Helix, Inc., Bozeman, MT, www.goldenhelix.com). The EMMAX-GRM analyses were conducted using an additive model that consisted of tests of association using single SNP (Kang et al., 2010) utilizing a genomic relationship matrix (VanRaden, 2008). The mixed model was specified as y = Xβ + Zμ + ϵ, where y is an n · 1 vector of the observed phenotypes, X is a n · q matrix of fixed effects, β is a q · 1 vector representing the levels of the fixed effects, and Z is a n · t matrix relating the instances of the random effect to the phenotypes. This assumes that and so that .In this study, Z is the identity matrix I, and K is the matrix of pairwise genetic relationships among samples (Kang et al., 2010). The variance components were estimated and the stratification among animals was accounted for and controlled using a GRM (VanRaden, 2008) computed from the genotypes. The individual EMMAX-GRM analysis for the Jersey population included state of origin (Oregon or Pennsylvania) as a covariate. The PNW and PNW + NE EMMAX - GRM analyses did not include any covariates as λ_GC = 1.0 and no other factors were significant in the prevalence of disease.

Statistical tests were conducted using an allele dosage model and 3 genotypic models (additive, dominance, recessive) in the PNW and PNW + NE populations with EMMAX (Kang et al., 2010). The allele dosage model was used to test for differences in allele frequencies between cases and controls. All 3 genotypic models were analyzed since the exact mode of inheritance at large effect QTL may not strictly be additive. The genotypic models were not applied to the Jersey population due to its small sample size. A QTL was defined by significant SNP(s) associated with Map tissue infection. Any SNP in linkage disequilibrium (LD) with an r² > 0.9 were considered to tag the same QTL. The regions surrounding significant QTL were further investigated to identify any putative positional candidate genes within 50 kb of the SNP(s) defining the most likely position of the QTL. An odds ratio (OR) for the major allele for each locus was used to determine the association between the locus and susceptibility to Map tissue infection (Morris and Gardner, 1988). The OR was calculated using the equation OR_M = (A · D)/(B · C) where OR_M denotes the OR for the major allele, A represents the number of major alleles in cases, D represents the number of minor alleles in controls, B represents the number of minor alleles in cases and C represents the number of major alleles in controls. Using this equation, an OR = 1 implies no association between the major allele for the locus and Map tissue infection. An OR > 1 implies the major allele is associated with increased risk, while an OR < 1 implies the major allele is associated with decreased risk of Map tissue infection.

The EMMAX models produced pseudo-heritability (percentage of phenotypic variance explained due to the fitted markers) estimates calculated as h² = σ²_g/(σ²_g + σ²_e) (Kang et al., 2010) for susceptibility to Map tissue infection in each population. The pseudo-heritability estimate along with the SE was reported for each population.

In the Jersey, PNW Holstein, and combined PNW and NE Holstein analyses performed with the BovineSNP50 genotypes, the Wellcome Trust Case Control Consortium (2007) threshold for GWAA was used where moderate evidence of association between a QTL and a trait was detected with P-values between 5.5 · 10^-5 and 5.5 · 10^-7 and strong evidence of an association was detected with P-values < 5.5 · 10^-7. In the combined PNW and NE Holstein meta-analysis where the genotypes were imputed to WGS level, the multiple testing correction method simpleM was used to determine the threshold for evidence of association (Gao et al., 2008). The simpleM method uses composite linkage disequilibrium to determine the correlation among the SNP and then infers the effective number of independent tests (M_eff) using principal components that explain 99.5% of the variation in linkage disequilibrium (Gao et al., 2008). The estimate of M_eff is then used in the Bonferroni correction to obtain the genome-wide significance threshold. For WGS with 16,063,340 SNP, we estimated the moderate significance threshold to be P < 1 · 10^-8.

RESULTS

Jersey Genome-wide Association Analysis

Seven QTL were identified as associated with Map tissue infection as determined by differences in allele frequencies between the Jersey cases and controls (Table 1). All QTL associated with Map infection were on BTA22 (Table 1). One QTL region on BTA22 encompasses the inositol 1, 4, 5-trisphosphate receptor type 1 (ITPR1) gene which was the most significant QTL identified in the Jersey population (P = 4.45 · 10^-6). Two other QTL were located within an intron of contactin-4 (CNTN4). All QTL identified in the Jersey cattle had OR < 1 (OR = 0.08 to 0.13). This indicates that the major allele for each SNP tagging a QTL (A nucleotide in rs41618810, rs41583533, rs41592357; C nucleotide in rs41640896, rs109685350; G nucleotide in rs41621239, rs110161711) is associated with a decreased susceptibility to Map tissue infection. The pseudo-heritability estimate for the Jersey population was 0.48 (SE = 1.20).

Pacific Northwest Holstein Genome-wide Association Analysis

Six QTL were associated with susceptibility to Map tissue infection in the PNW Holstein cattle GWAA (Table 2). The dominance model identified QTL on BTA8 and BTA14 and 2 on BTA21, while the recessive model identified a single QTL on BTA3. The QTL on BTA3 was located within the breast cancer anti-estrogen resistance 3 (BCAR3) gene. A second QTL on BTA8 and a QTL on BTA14 were associated with susceptibility to Map tissue infection using both the allelic and dominance models. The most significant (P = 8.6 · 10^-8) QTL identified in the PNW population was on BTA21 and was identified using the allelic, additive, and dominance models. There were no positional candidate genes located in the 100 kb regions surrounding the lead-SNP for these 3 QTLs. Four QTL (located on BTA3, BTA8, and BTA14) identified in the PNW Holstein GWAA had OR < 1 (OR = 0.31 to 0.48) indicating that the G allele on BTA3 (rs110936943), the G alleles on BTA8 (rs29010355 and rs42454958), and the T allele on BTA14 (rs29016507) decreased the risk for becoming Map tissue infected. Two QTL (both located on BTA21) had OR > 1 (OR = 2.53 to 3.54), indicating that the animals with T (rs109376203) and G (rs41585301) alleles were 2 to 3 times more likely to be Map tissue infected. The pseudo-heritability estimate for the PNW Holstein population was 0.12 (SE = 0.25).

Combined Pacific Northwest and Northeast Holstein Genome-wide Association Analyses

In the combined PNW and NE Holstein population meta-analysis using the BovineSNP50 data, 6 QTL were associated with Map tissue infection (Table 3). The allelic model identified 2 QTL (BTA10 and BTA12), the dominance models identified 3 QTL (BTA8, BTA14, and BTA16), and the recessive model identified 1 QTL (BTA16). The QTL on BTA8 and BTA14 had previously been identified to be associated with Map tissue infection in the PNW population (Table 2). Two positional candidate genes, feline leukemia virus subgroup c cellular receptor family member 2 (FLVCR2) and chromosome 10 open reading frame 1 (C10H14ORF1), were identified within the region harboring the QTL on BTA10. RAS p21 protein activator 3 (RASA3) and an uncharacterized protein (MGC134473) were identified in the 100 kb region surrounding the QTL on BTA12, which was the most significant QTL detected in the combined PNW and NE Holstein population (P = 2.33 · 10^-5). The QTL region on BTA16 contained 3 positional candidate genes: microtubule affinity regulating kinase 1 (MARK1), mitochondrial amidoxime-reducing component 2 (MARC2), and chromosome 1 open reading frame 115 (C16H1orf115). One QTL on BTA14 was found by both the allelic and dominance models and another BTA16 QTL, containing 3 associated SNP, was found by both the additive and dominance models within introns of CDC42 binding protein kinase ɑ (CDC42BPA). Two of the 3 SNP (rs42396640 and rs42397377) had been previously associated with Map tissue infection in the NE Holsteins by Settles et al. (2009). Five QTL identified in the PNW and NE BovineSNP50 meta-analysis had OR < 1, indicating the major alleles (G nucleotide in rs42454958; C nucleotide in rs29014131; T nucleotide in rs29016507, rs110520204; A nucleotide in rs42396640) were associated with decreased susceptibility to Map tissue infection (OR = 0.42 to 0.63). While 1 QTL (tagged by rs43170466) on BTA10 had an OR > 1, indicating that the major allele (C) was associated with increased susceptibility to Map infection (OR = 1.86). The pseudo-heritability estimate for the combined PNW and NE Holstein population was 0.07 (SE = 0.14).

In the combined PNW and NE Holstein population meta-analysis using genotypes imputed to WGS level, none of the QTL met the significance threshold of P < 1 · 10^-8. The QTL closest to the significance threshold was on BTA1 (rs381822862; P = 9.34 · 10^-8) and was located near the butyrylcholinesterase (BCHE) gene. The pseudo-heritability estimate for the combined PNW and NE Holstein population meta-analysis using genotypes imputed to WGS level was 0.05 (SE = 0.09), which was slightly less than the pseudo-heritability estimate obtained for the BovineSNP50 data.

DISCUSSION

Three populations (Jersey, PNW Holstein and combined PNW and NE Holstein) were investigated for SNP associations with susceptibility to Map tissue infection using 4 (allelic, additive, dominance, and recessive) GWAA models in the PNW and PNW + NE populations, and the allelic dosage model in the Jersey population. These populations were investigated to identify new loci associated with Map tissue infection and to validate associations that had been identified in 8 previous cattle studies investigating loci associated with bovine paratuberculosis (Settles et al., 2009; Minozzi et al., 2010; Pant et al., 2010; Kirkpatrick et al., 2011; van Hulzen et al., 2012; Minozzi et al., 2012; Alpay et al., 2014; Zare et al., 2014). This included a meta-analysis by Minozzi et al. (2012) which included cattle from both the Settles et al. (2009) and the Minozzi et al. (2010) analyses.

Only 3 QTL associated with susceptibility to bovine paratuberculosis have been found in more than 1 study. One QTL consisting of 2 SNP on BTA1, rs29012843 [(3 megabases (Mb)] and rs29012842 (3 Mb), was identified in a GWAA by Settles et al. (2009; P = 3.26 · 10^-6) and were also found by the meta-analysis performed by Minozzi et al. (2012; P = 4.82 · 10^-6) although it is important to note that both studies shared some of the same animals. There were no positional candidate genes located within 50 kb of this QTL. On BTA12, 2 QTL, rs42743330 (70 Mb) and rs41665666 (72 Mb), were identified by Minozzi et al. (2010; P = 3.55 · 10-6; P = 1.44 · 10^-6) and were found by the meta-analysis of Minozzi et al. (2012; P = 2.66 · 10-5; P = 2.88 · 10^-5) which also included cattle that were shared between studies. Neither QTL had positional candidate genes located within 50 kb of the associated lead-SNP. While these 3 QTL were not validated in independent populations, they are the only QTL to have been identified in multiple analyses, to date.

Two SNP (rs42396440 and rs42397377) identified by Settles et al. (2009; P = 2.57 · 10^-5) were identified by the combined PNW and NE Holstein meta-analysis in this study. While these markers were not validated in an independent population, the identification of these markers in the combined PNW and NE Holstein meta-analysis suggests that these previously identified markers associated with Map tissue infection warrant further investigation. No other QTL found in the current analyses lie within 50 kb or in LD (all r² < 0.9 threshold) with the lead-SNP that tagged the QTL identified in these 8 previous studies. A QTL on BTA14 (rs29016507), identified in both the PNW Holstein and the combined PNW and NE Holstein meta-analysis, is located approximately 1.4 Mb away from 2 SNP (ss61521480, P = 7.8 · 10^-9; rs42413954, P = 6.8 · 10^-8) that were previously identified by Pant et al. (2010). Another QTL on BTA8 (rs42454958), also identified in both the PNW Holstein and the combined PNW and NE Holstein meta-analysis, is located approximately 1.4 Mb away from 1 SNP identified by Minozzi et al. (2010; rs43161947, P = 7.02 · 10^-5).

Comparisons across studies were hampered by the different methods used to diagnose cattle with Johne's disease, the definition of control cattle, and the version of the Illumina Bovine SNP50 BeadChip used. Unfortunately, later versions of the BovineSNP50 BeadChip resulted in the elimination of SNP that had previously been shown to be associated with Map tissue infections in early studies, making it more difficult for those initial studies to be directly validated. In Settles et al. (2009), the “tissue” definition of the Johne's disease phenotype was most consistent with the Johne's disease phenotype used in the PNW and Jersey cattle of this study. Settles et al. (2009) cultured tissues taken from ileocecal nodes to determine if an individual was Map infected, whereas, the current study conducted qPCR from the ileocecal lymph nodes to determine if an individual was Map infected as an indicator of the presence of Johne's disease. Others (Pant et al., 2010; Kirkpatrick et al., 2011) have used ELISA (blood or milk) or fecal culture results to indicate the presence of Johne's disease.

The specificity of diagnostic tests across sample sources and methods are similar (Whitlock et al., 2000; McKenna et al., 2005; White et al., 2016). Reports of the sensitivities of diagnostic test for ELISA, cultures and qPCR of serum, milk, feces, and tissues are difficult to compare as methodologies and samples tested differ. The test sensitivity of ELISA tests based on fecal culture positive cows are 43% to 94% (Colgrove et al., 1989, Collins et al., 1991; Sockett et al., 1992). However, when sensitivities are based on the level of fecal shedding (as cattle shedding a higher number of organisms are more likely to be positive by ELISA), sensitivity drops to 15% for low Map shedders and 88% for cattle shedding high amounts of Map (Dargatz et al., 2001). The sensitivities of ELISA and fecal culture when compared with tissue culture positive cattle, drop to 6.9% to 16.9% for ELISA and 13.9% to 27.8% for fecal culture (McKenna et al., 2005). The testing sensitivity for qPCR ranges from 87% to 95% although the source of the sample affects that sensitivity as only 7% of Map tissue positive samples were identified as positive from fecal samples (Whitlock et al., 2000; Erume et al., 2001; White et al., 2016).

The lack of detection of common QTL across the studies could easily be due to differences in sensitivity of the various diagnostic methods and the sources of samples used to detect Map infection in its early stages in this and the previous 8 studies investigating loci associated with bovine paratuberculosis. The selection of a standard diagnostic method and phenotype for bovine paratuberculosis would make identifying QTL associated with susceptibility to Map infection easier. Sensitivity of detection of Map would favor the use of diagnostic testing of tissues for research studies, but this method of detection is costly, is only practical postmortem, can be taken at a single time point, and is difficult to do on large numbers of animals. The use of less sensitive diagnostic methods facilitates studies with larger numbers of animals over a time continuum, but is more likely to have control animals that are actually infected with Map.

Unfortunately, the choice of diagnostic method and sample source results in trade-offs in the power of the analyses. More sensitive diagnostics and samples from tissues increase the power of the study by reducing false negative controls, but also results in a reduction of power because of subsequent small sample sizes. The use of less sensitive diagnostic methods and samples that are more convenient is conducive to larger sample sizes, thereby increasing power, but this power may be lost or partially lost due to inaccurate diagnosis of control samples. A standardized phenotype is therefore dependent on the diagnostics and the samples taken; until a highly sensitive non-invasive diagnostic test is developed that identifies cattle with early stage disease, the struggle of identifying and validating QTL associated with bovine paratuberculosis will continue to be hampered.

Many of the QTL associated with susceptibility to Map tissue infection found in this study were located near positional candidate genes that have functions related to signal transduction pathways and the immune system. Given the nature of Map infection, several of these positional candidate genes (CDC42BPA, BCAR3, and RASA3) also have functions that could be related to how Map enters the host cells or how the immune responses are initiated due to Map infection. Once Map is ingested, there are 2 main mechanisms that the organism uses to cross the intestinal epithelium: endocytosis by microfold cells (M cells) or uptake by enterocytes using the cell division cycle 42 (CDC42) - ras homolog family member A (RhoA) pathway (Bannantine and Bermudez, 2013). Many positional candidate genes are associated with the regulation of GTPases, like CDC42 and RhoA, including BCAR3, CDC42BPA, and RASA3. The BCAR3 gene is involved in signal transduction pathways and has been shown to help activate small GTPases (Rufanova et al., 2009). CDC42 binding protein kinase ɑ functions as a downstream mediator of the GTPase CDC42 and is thought to contribute to transferrin receptor-mediated iron uptake (Cmejla et al., 2006). The RASA3 gene is a GTPase that is involved with the Ras signaling pathway that regulates cell proliferation, growth, migration and survival, as well as the regulation of other GTPases (Blanc et al., 2012). Modification of the function of GTPases, such as CDC42 or RhoA, could impact host susceptibility by altering the pathway(s) that Map uses to cross the intestinal epithelium and infect the host. For example, a study in mice found that when CDC42 was deleted in a Cre transgenic mouse model, mouse pups died at 3 wk of age due to infections in the respiratory tract resulting from a loss in the hosts' ability to control infections (Lee et al., 2013).

Two positional candidate genes are involved with signal transduction pathways that regulate cell migration: BCAR3, and ITPR1. Cellular migration plays an important role throughout an organism's life through processes such as immune surveillance and wound repair (Franz et al., 2002). Without the ability to move, immune cells such as macrophages would be unable to infiltrate infected tissues or lymph nodes to stimulate an immune response. Given that Map can disseminate throughout the host after the initial infection, mutations to cell migration could also impact the ability of infected cells to spread through the body, as well as the ability of additional immune cells to travel to the site of infection to further the immune response (Bannantine and Bermudez, 2013).

Two candidate genes are involved with regulation of the WNT signaling pathway: ITPR1 and MARK1. The WNT signaling pathway is activated during fibroblast proliferation that occurs in tissues being repaired after damage caused by inflammatory lesions (Cheon et al., 2004). A distinguishing feature of Map infection is the formation of granulomas which can form as the host's macrophages attempt to contain the initial infection and as an abnormal proliferation or hyperplasia of fibroblasts (Ackermann, 2007). Any mutation that potentially perturbs the WNT pathway could hinder or enhance a host's immune response to inflammatory lesions found in the intestines of infected cattle.

Two candidate genes (ITPR1 and RASA3) regulate platelets, which function in the innate and adaptive immune responses through secretion of proinflammatory cytokines (Conglei et al., 2012). Inositol 1, 4, 5-trisphosphate receptor type 1 regulates platelet aggregation, which involves platelet-to-platelet adhesion after an initial response of platelets to a site of injury, and platelet activation which occurs after aggregation (Rumbaut and Thiagarajan, 2010). Conversely, RASA3 is a known inhibitor of platelet activation (Stefanini et al., 2015). Platelets become activated when several different agonists and collagens bind to cell surface receptors to ultimately increase the concentration of calcium in the cytoplasm of the platelets and alter their morphology. Platelets can also be activated by the binding of bacteria which initiates secretion of cytokines as well as antimicrobial peptides (Cox et al., 2011). Mutations that alter the function or regulation of platelets could hinder the ability of the host to fight a Map infection by limiting the inflammatory response or the recruitment of lymphocytes to infected areas.

Only 1 (of 6) QTL associated with susceptibility to Map tissue infection in the PNW Holsteins was located near a positional candidate gene. The lack of positional candidate genes located near QTL, however, does not suggest that the QTL play no role in altering gene expression. The majority of causal mutations lie within cis-regulatory regions of the genome (3' and 5' untranslated regions, promoters, enhancers, and introns) that are located some distance (50 kb to over 1 Mb) away from the transcriptional start site of a gene (Koufariotis et al., 2014). It is, therefore, expected that the majority of DNA variants that will result in gene expression differences will be found in non-coding DNA regions that affect distant genes rather than in coding regions.

The positional candidate genes identified in this study were compared to positional candidate genes that were previously identified in 20 Crohn's disease GWAA to determine if there was any overlap. Crohn's disease is a chronic inflammatory disease of the gastrointestinal tract in humans that shares many of the same symptoms as Johne's disease (Hawkey et al., 2012). However, Crohn's is a multifactorial disease, as there are several pathogens and/or physiological stressors that trigger the disease. There is concern that Map has the potential to spread from infected livestock into the human population. Given that Map has been isolated from both healthy individuals and Crohn's patients, it might be an opportunistic bacteria living within the human gastrointestinal tract (Sechi et al., 2005). None of the positional candidate genes identified in this study have previously been associated with risk of Crohn's disease. However, 2 genes associated with risk of Crohn's disease, breast cancer anti-estrogen resistance 1 (BCAR1) and inositol 1, 4, 5-trisphosphate receptor type 3 (ITPR3) were closely related to positional candidate genes that were identified in this study (BCAR3 and ITPR1). Many genes previously linked to susceptibility to Crohn's disease have also been linked to other mycobacterial infections in humans (Behr and Divangahi, 2015; Nabatov, 2015). None of the positional candidate genes identified in this study were associated with other mycobacterial diseases. However, identifying genes associated with multiple mycobacterial infections would provide researchers with the ability to further understand the differences in pathogenesis of each mycobacterium as well as identify populations at-risk for developing mycobacterial diseases.

Heritability estimates for susceptibility to bovine paratuberculosis of 0.06 to 0.16 have been reported (Koets et al., 2000; Gonda et al., 2006). Several factors contribute to the variation in reported heritability estimates including level of Map exposure, the accuracy of the phenotype as a result of the diagnostic test and samples used, herd management, sample size, and contemporary groups. As was previously discussed, variations in the accuracy of the phenotypes across studies is likely, and other factors such as sample size have also contributed to differences in heritability estimates. In the PNW Holstein and the two combined PNW and NE Holstein population analyses, the pseudo-heritability was relatively low; this could be due to increased environmental variation between the cattle populations, possible age differences in the cattle since the age of the PNW Holstein cows was unknown, and potential differences in Map exposure. The large pseudo-heritability in the Jersey cattle could be due to the high level of Map exposure for the 9 steers or it could be an imprecise estimate as the SE is quite large. The slightly lower pseudo-heritability for the PNW and NE Holstein meta-analysis using imputation data seems unintuitive as traditionally, having more markers available increases the heritability estimate (Visscher et al., 2008). This lower pseudo-heritability may be due to errors during the imputation process, which was performed using a reference population of multiple breeds that might have contributed to the error. Effective use of genetic selection is still possible for traits with low heritability estimates. In fact, a recent study on the impact of genetic section on US Holstein cattle found that lowly heritable traits such as daughter pregnancy rate and somatic cell score had greater response to genomic selection compared to moderately heritable traits like protein yield (García-Ruiz et al., 2016).

The population sizes used in this study were smaller than some of the previously published studies that have identified QTL associated with bovine paratuberculosis. However, the current study is unique in the definition of cases and controls through the use of qPCR detection of Map in tissue samples rather than by ELISA or culture of serum, milk, feces, or a combination of these. Tissue culturing/testing is considered to be the gold standard in diagnostic testing due to its increased sensitivity and ability to detect animals in earlier stages of the disease that can be missed with other diagnostic methods and decreases the likelihood of false negative controls (Pavlik et al., 2000). Additionally, only 1 other GWAA investigating susceptibility to Map infection in Jersey cattle has been published (Zare et al., 2014), so this study provides additional information about susceptibility to Map infection in that breed.

The objectives of this study were to validate previously identified loci within Jersey and Holstein cattle, as well as to identify new loci associated with susceptibility to Map tissue infection. The QTL on BTA16 previously identified in the Settles et al. (2009) NE Holstein GWAA, was also identified by the combined PNW + NE Holstein meta-analysis. Although these 2 GWAA did share some animals in common, the combined analysis did not result in a GWAA with a smaller P-value than was initially reported by Settles et al. (2009). This may have been due to a lack of power due to small sample sizes, differences in the genetics underlying the 3 populations, a weak association on BTA16 with Map tissue infection or possibly a false positive. Sixteen new QTL were found in the Jersey and Holstein populations: 7 loci were identified in Jersey cattle and 9 loci were identified in Holstein cattle. The QTL on BTA8 and on BTA14 that were identified in the individual PNW population and the combined PNW + NE population, although not replicated in independent populations, do suggest that these QTL warrant further investigation.

The identification and validation of QTL associated with susceptibility to Map tissue infection in cattle could be advanced through the use of a consistent diagnostic method and phenotype classification. Future studies need to be performed to determine if these newly identified QTL are affecting the positional candidate genes near-by or other genes through cis- or trans- regulatory effects. If the loci are located within genes, experiments such as CRISPR knockouts or luciferase assays could be performed to determine how the mutation is affecting gene expression. For loci located in intergenic regions of the genome, chromosome conformation capture approaches could be used to determine if the loci are located in regulatory elements and which genes these regulatory elements are interacting with. Once QTL are confirmed to be associated with Map infection, their addition to commercially available genotyping assays would allow producers to select cattle that are less susceptible to the disease, ultimately reducing the spread of the disease and preventing further economic losses.

LITERATURE CITED