Single Nucleotide Polymorphisms in the Bovine TLR2 Extracellular Domain Contribute to Breed and Species-Specific Innate Immune Functionality.

Recent evidence suggests that several cattle breeds may be more resistant to infection with the zoonotic pathogen Mycobacterium bovis. Our data presented here suggests that the response to mycobacterial antigens varies in macrophages generated from Brown Swiss (BS) and Holstein Friesian (HF) cattle, two breeds belonging to the Bos taurus family. Whole genome sequencing of the Brown Swiss genome identified several potential candidate genes, in particular Toll-like Receptor-2 (TLR2), a pattern recognition receptor (PRR) that has previously been described to be involved in mycobacterial recognition. Further investigation revealed single nucleotide polymorphisms (SNP) in TLR2 that were identified between DNA isolated from cells of BS and HF cows. Interestingly, one specific SNP, H326Q, showed a different genotype frequency in two cattle subspecies, Bos (B.) taurus and Bos indicus. Cloning of the TLR2 gene and subsequent gene-reporter and chemokine assays revealed that this SNP, present in BS and Bos indicus breeds, resulted in a significantly higher response to mycobacterial antigens as well as tri-acylated lipopeptide ligands in general. Comparing wild-type and H326Q containing TLR2 responses, wild-type bovine TLR2 response showed clear, diminished mycobacterial antigen responses compared to human TLR2, however bovine TLR2 responses containing H326Q were found to be partially recovered compared to human TLR2. The creation of human:bovine TLR2 chimeras increased the response to mycobacterial antigens compared to the full-length bovine TLR2, but significantly reduced the response compared to the full-length human TLR2. Thus, our data, not only present evidence that TLR2 is a major PRR in the mammalian species-specific response to mycobacterial antigens, but furthermore, that there are clear differences between the response seen in different cattle breeds, which may contribute to their enhanced or reduced susceptibility to mycobacterial infection.

Introduction

Early recognition of pathogens and activation of the innate immune response is critical in determining the outcome of infection. Pathogens are recognized by evolutionary conserved, germline-encoded and cell-surface expressed pattern recognition receptors (PRR) (1). One of the major groups of PRRs are the toll-like receptors (TLR), ten of which have been discovered common to both bovine and human immune systems (2–5), and are expressed on myeloid antigen presenting cells such as macrophages (MØ). Each of these TLR has evolved to recognise specific pathogen associated molecular patterns (PAMP) – and docking of these ligands to the extracellular domain (ECD) of the TLR triggers an intracellular signalling cascade (6). This signalling is mediated by the recruitment of adaptor proteins to the intracellular Toll/interleukin-1 receptor homology (TIR) domain of the TLR (7, 8). While the amino acid sequence of the ECD varies considerably between species (9), the TIR domain is highly conserved (7, 9).

Although TLR4, and possibly TLR9 and TLR8 are capable of sensing mycobacteria (10–12), recognition of mycobacterial antigens by TLR2 has been shown to be a key factor in determining disease progression in tuberculosis (12–15) (14, 16–18). TLR2 recognizes a range of mycobacterial cell-wall antigens including the 19 kDA mycobacterial lipoprotein and glycolipids such as lipoarabinomannan (LAM), lipomannan (LM), the 38 kDA antigen, LprG lipoprotein and phosphatidylinositol mannoside (PIM) (11, 19, 20). Additionally, secreted mycobacterial proteins have been shown to act as TLR2 ligands. For instance, the well-studied virulence factor ESAT-6 has been shown to inhibit NF-κB activation in a TLR2-mediated manner (21). Similarly, recombinant antigen TB10.4 from M. bovis, encoded by a subfamily of ESAT-6 and also present in M. tuberculosis and M. bovis BCG (22), has been shown to induce pro-inflammatory cytokine production via a TLR2 mediated NF-κB pathway (23).

Engagement of TLR2 by mycobacteria can result in complex outcomes that can either promote or inhibit a protective inflammatory response. For example, TLR2 activation can result in the upregulation of pro-inflammatory cytokines including IL-1β, TNFα and IL-6 as well as reactive oxygen species (ROS) production through activation of NF-κB and MAPK pathways (24–26). Furthermore, TLR2 signalling regulates the expression of iNOS, which interacts with NADPH oxidase and ROS to create reactive nitrogen intermediates (27, 28) and can induce autophagy (29) - all essential mechanisms in cellular host defence to mycobacterial infection (30, 31). On the other hand, engagement of TLR2 by the secreted mycobacterial antigens heat-shock protein (hsp) 60 and proline-proline-glutamic acid (PPE)18 initiates an anti-inflammatory response via IL-10 production (32, 33). In addition, ESAT-6 has been shown to inhibit IL-12 production upon binding to TLR2 in RAW cells (21). Further immunomodulatory mechanisms of some mycobacterial species include the inhibition of MHC II receptor expression and therefore presentation of mycobacterial lipoprotein antigens (34, 35). Thus, mycobacteria can utilise TLR2 to manipulate the immune outcome and promote survival (15, 36) (13, 14, 19, 21, 29, 36).

Given the importance of TLR2, it is perhaps not surprising that polymorphisms within the TLR2 gene have been shown to affect the immune response to mycobacterial infection in both humans and cattle (18, 29, 32, 33, 37, 38) (12, 14),. In humans, the R753Q and R677W mutations have been shown to influence susceptibility to tuberculosis in case-control studies (37, 39). Furthermore, the R753Q SNP was also shown by Pattabiraman et al. (40) to result in impaired NF-κB signalling and decreased cytokine responses to M. smegmatis in a murine and Human Embryonic Kidney (HEK) 293 cell assay. Similar as in the human system, several SNPs have also been identified within the bovine TLR2 orthologue. Some of these have been indicated to contribute to disease susceptibility in cases of paratuberculosis (41), clinical mastitis (42, 43), and bovine tuberculosis (44–46). SNPs occurring between cattle breeds may reflect natural variation arising from host-pathogen co-evolution in different geographical environments, but may also be influenced by intensive selective breeding that favours meat or milk production (9, 47). Characterisation of SNPs arising in TLR2 among commercially important breeds is therefore a crucial contribution towards understanding immune fitness variation among herds and potentially, disease sensitivity. While some SNPs within the bovine tlr2 gene have previously been identified, studies of their functional relevance are often lacking. In the present work, we compare the tlr2 sequences of two breeds of global importance to the farming/dairy industry, Brown Swiss (BS) and Holstein Friesian (HF) and determine the genotypic frequency of a H326Q SNP that occurs in the ligand-binding region of the ECD. We report on the species- and breed- specific phenotypic consequence of this mutation in response to canonical TLR2 ligands and in the context of challenge with M. tuberculosis and M. bovis. We present evidence that TLR2 is a major PRR in the mammalian species-specific response to mycobacterial antigens, that there are clear differences in species-specific responses that may influence disease susceptibility and host tropism of Mycobacterium tuberculosis complex members.

Methods

Mycobacterial Strains and Culture

The following strains were used for all assays as indicated: M. bovis strain AF2122/97 and M. bovis BCG Tokyo; M. tuberculosis H37Rv. Strains were cultured in BD Middlebrook 7H9 media (BD Difco™, USA), supplemented with OADC and 0.05% Tween 80 (Sigma-Aldrich, UK). For M. bovis culture only, growth media was additionally supplemented with sodium pyruvate; for M. tuberculosis culture only, media was additionally supplemented with glycerol (Sigma-Aldrich, UK).

Brown Swiss Genome Assembly and QC

Four lanes of 10x Genomics data were generated using an Illumina HiSeq. This provided ~1200 million reads with a 150bp read length. The data was assembled de novo using the 10x Genomics assembly tool Supernova v2.0.1, using the option to output a phased genome, which generates two homologous assemblies where each scaffold in the assembly has one of the paternal or maternal haplotypes (48). Using one of the haplotypes, we used RepeatMasker, using ‘bovidae’ as the repeat library to soft-mask repeat regions in the assembly (49). We tested the assembly for completeness using several quality control metrics. First, BUSCO (v5.0.0) was used to identify the number of mammalian single-copy orthologs present in the genome (using the database mammalia_odb10) (50). Using the Kmer Analysis Toolkit, we examined the content and distribution of kmers in the assembly to that of the reads using the tool ‘kat comp’, ignoring the first 16 bases of each R1 read in order to omit 10x barcodes from the analyses. We then visualised the distribution of kmers using the tool ‘kat plot’ (51).

Brown Swiss Genome Annotation and Gene Expansion

The LiftOff tool was used to lift over the B. taurus genome annotation (ARS-UCD1.2.103) to that of the Brown Swiss, noting any genes that were unmapped in Brown Swiss (52). We then identified genes in the Brown Swiss genome annotation that were present in multiple copies (compared to B. taurus ARS-UCD1.2.103) and downloaded the functional information for these genes using ENSEMBL BioMart database (53). For genes annotated as being duplicated three or more times and where a gene name or description was missing, we downloaded the coding sequence for that gene and used megablast to search the NCBI nr database and noted the highest scoring hit (with 100% query coverage and an E-value equal to 0), that contained functional information about the gene (i.e., disregarding hits for clones, BACs, isolates, etc.).

Study Population for Candidate Gene TLR2 Analysis

Complementary DNA (cDNA) from TLR2 mRNA sequences isolated from PBMCs from six clinically healthy pedigree BS and four HF were investigated for sequence comparison. All cows were female and housed under Home Office License PPL7009059. Age and lactation cycle of sampled cows were recorded. After identification and selection of candidate SNP H326Q within the coding sequence (CDS) of TLR2 in the BS samples, further cDNA was generated from nine BS and 13 HF PBMC samples from the above-mentioned farms to be assessed for the presence of this SNP. Additionally, 17 Sahiwal and 17 Boran cDNA samples were kindly provided by Drs Thomas Tzelos and Tim Connelley (both The Roslin Institute, University of Edinburgh, UK) for TLR2 SNP analysis.

PBMC Isolation and Macrophage Maturation

Blood for peripheral blood mononuclear cells (PBMC) isolation and subsequent macrophage (MØ) generation was collected from clinically healthy pedigree HF and BS cows. All procedures were carried out under the Home Office license (PPL7009059) approved by the RVC’s AWERB Committee. Blood was drawn into sterile glass vacuum bottles containing 10% v/v acid citrate and centrifuged to separate the buffy coat, which was then diluted in phosphate buffered saline (PBS, Thermo Fisher Scientific, UK) with 2% foetal calf serum (FCS, Sigma-Aldrich, UK). Lymphoprep (d = 1.077 g mL-1, StemCell Technologies ™, Canada) was underlaid and PBMC were isolated by density gradient centrifugation. PBMC were washed and red blood cells lysed, resultant cells were cultured in RPMI 1640 with GlutaMAX™ (Gibco, UK), supplemented with FCS (Sigma-Aldrich, UK) and penicillin/streptomycin (Thermo Fisher Scientific, UK) at 37 °C with 5% CO2. To derive MØ, PBMC were incubated additionally with 10% filtered L929 supernatant. Media was replaced after three days and adherent cells were harvested after 6 days. Matured MØ were phenotypically assessed by flow cytometry using antibodies specific for bovine cell surface expressed markers including CD14, CD16, CD32, CD11b, CD80, CD163 and MHC class II. All antibodies and isotype controls were purchased from Bio-Rad, UK ( ). Data were acquired using a FACS Calibur (BD Biosciences, UK) with Cell Quest Pro software (BD Biosciences, UK), counting 10,000 events. Data was exported as FSC files and analysed using FlowJo V10 (FlowJo LLC, USA).

Table 1

Antibodies used in this study.

Antibody	Clone	Label	Dilution
Mouse IgG1 negative control	–	FITC	1:20
Mouse IgG1 negative control	–	RPE	1:20
Mouse IgG2a negative control	–	FITC	1:20
Mouse IgG2b negative control	–	RPE	1:20
Mouse IgG2a anti Human CD 14	TÜK4	RPE	–
Mouse IgG2a anti Human CD 16	KD1	FITC	–
Mouse IgG2b anti Bovine CD 11b	CC126	FITC	–
Mouse IgG1 anti Bovine CD 32	CCG36	RPE	–
Mouse IgG1 anti Bovine CD 80	IL-A159	FITC	–
Mouse IgG1 anti Human CD 163	EDHu-1	RPE	–
Mouse IgG1 anti Bovine MHC II: DR	CC108	RPE	1:25
Human IgG anti Bovine CD282		AF647
Mouse IgG anti Human CD282		AF647
HuCAL Fab-dHLX-MH		AF647

RNA Extraction and cDNA Synthesis

For sequencing and cloning of the tlr2 gene, total RNA was isolated from bovine PBMC using a RNeasy Mini kit (Qiagen, UK) with on-column gDNA digestion (DNAse from Ambion, UK) according to the manufacturer’s instructions. Total RNA was reverse transcribed to cDNA using iScript cDNA Synthesis Kit (Bio-Rad, UK) according to the manufacturer’s instructions. cDNA samples were checked for concentration and purity via the A260/280 ratio using a Nanodrop ND-1000 Spectrophotometer (Thermo Fisher Scientific, UK). A RT-negative control of each RNA sample was analysed by polymerase chain reaction (PCR) for expression of β-actin to confirm lack of gDNA contamination of RNA.

Cloning and Sequencing of tlr2

Full length tlr2 sequences were amplified using primers designed using Primer-BLAST (NCBI, USA) and Primer3Plus ( ) (187) based on bovine genome Bos_taurus_UMD_3.1.1 (NCBI, USA) as reference. PCR reactions were performed using Easy-A High-Fidelity PCR master mix (Agilent Technologies, UK) on a Mastercycler Pro Thermal Cycler (Eppendorf, UK). Following amplification, the PCR products were purified using the MinElute PCR purification kit (Qiagen, UK) according to the manufacturer’s protocol. For TA cloning of full length bovine tlr2 into the pDrive plasmid (Qiagen, UK) or pGEM-T plasmid (Promega, UK), manufacturer’s instructions were followed. Plasmids were transformed into NEB® 10-beta competent E. coli (New England Biolabs, UK) and colony PCR performed to verify the presence of the tlr2 insert. Sequences were validated using Sanger sequencing; DNA sequencing was performed by DNA Sequencing & Services (MRC I PPU, School of Life Sciences, University of Dundee, Scotland, www.dnaseq.co.uk) using the corresponding vector-specific sequencing primers. A minimum of 4 different plasmids derived from unique colonies were sequenced for each animal. Bovine TLR2 contigs were assembled using the bovine TLR2 sequence Bos taurus Hereford breed NM_174197.2 (NCBI, USA) as a reference. TLR2 nucleotide sequences were translated into corresponding amino acid sequences and aligned in CLC Main Workbench V7.6.4 (CLCbio, Denmark). Detected SNPs were classified according to the dbSNP database (NCBI, USA). Analysis of SNP data was curated and as at present, the NCBI dbSNP database contains only human data, allele changes, CDS position and function of each SNP was manually cross-checked and referenced against Ensemble and the European Variation Archive (EVA).

Table 2

Primer names for PCR reactions.

Primer name	Sequence	GC (%)	Tm (°C)
bovine TLR2 cds forward	CACCATGCCGCGTGCTTTGTGGACA	60	67.9
bovine TLR2 cds reverse	CGAAGGGTCCTAGGACCTTATTGCAGCTCTC	54.8	70.8
bovine TLR2 int forward	CAGTGCTCAAAATCTGCAGATA	40.9	56.5
bovine TLR2 int reverse	GTACTCATTTCACTGATGGATGC	43.5	58.9

All primers were designed using Primer-BLAST (NCBI) and Primer3Plus.

Expression of TLR2 by HEK 293T (SEAP) Cells

The Human Embryonic Kidney (HEK) 293T cell line, stably transfected with the gene encoding NF-κB inducible secreted embryonic alkaline phosphatase (SEAP), was cultured SEAP media as described recently (54). Bovine or human TLR2 sequences of interest were cloned into pcDNA3.3 TOPO TA mammalian expression vector (Thermo Fisher Scientific, UK), and transfected into SEAP HEK cells using TurboFect transfection reagent (Thermo Fisher Scientific, UK) according to the manufacturer’s protocol. After 24 h, cells were counted and re-seeded in triplicates in a 96-well plate to rest for another 24 h before either evaluation of transfection efficiency or stimulation with various agonists. Surface expression of TLR2 was confirmed by flow cytometry using either human anti-bovine CD282 (anti-TLR2) antibody or mouse anti-human CD282 antibody, both conjugated with Alexa Flour 647 (Bio-Rad, USA). HuCAL Fab-dHLX-MH Ab (Bio-Rad, USA) was used as a negative control as recommended by the manufacturer.

Stimulation Assays

We recently described the SEAP assay and the correlation between SEAP and CXCL8 production in detail (54), see also . In brief, SEAP HEK cells expressing TLR2 constructs, or PBMC-derived bovine MØ were specifically stimulated using 100 ng mL-1 diacylated FSL-1 (a TLR2/TLR6 antagonist; In vivogen, USA); 1 mg mL-1 of Pam3CSK4 (a TLR2/TLR1 antagonist, In vivogen); 19 kDa lipoprotein antigen (represents Rv3763 or Mb3789) (EMC microcollections, Germany), M. bovis BCG (MOI of 10) or MTB H37Rv at (MOI of 5). The direct NF-κB activator phorbol 12-myristate 13-acetate (PMA, Sigma-Aldrich, UK) was used at 200 ng mL-1 as a positive control. Sterile culture media and mock transfection with an empty vector plasmid were used as negative controls. Supernatants were harvested for CXCL8 ELISA as described by Cronin et al. (55) for bovine MØ, or Quantikine Human CXCL8 ELISA kit (R&D Systems, USA) for human MØ, or colorimetric SEAP assay from HEK cells. SEAP activity as induced by NF-kB activation was assessed by the addition of pre-warmed Quantiblue reagent (In vivogen, USA) followed by incubation for 24 h at 37°C. For both assays, optical density was measured using a SpectraMax M2 spectrometer (Molecular Devices, UK).

Results

We used the 10x Genomics tool Supernova to assemble a haplotype-phased Brown Swiss cow genome. After trimming, the mean read length was 138.5 bases, with a raw read coverage of 59.64. The genome size was calculated to be 2.66Gb (of which 98.08% of the genome was present in scaffolds >10Kb) with a contig N50 of 523Kb and a scaffold N50 of 26Mb ( ). We used BUSCO to identify reconstructed single copy orthologs. From 9226 BUSCO orthologs 95.2% were recovered, 93.1% were single-copy, 2.1% were duplicated, 1.5% were fragmented, and 3.3% were missing. We used KAT to examine the distribution of kmers in the assembly to that in the reads ( , red distribution). The genome shows a normal distribution of kmers found once in the assembly and once in the reads. It also shows a number of low frequencies kmers (black peak between 0-10) present in the reads, but not in the assembly that most likely represent sequencing errors removed from the assembly. This continues into a shallow distribution of kmers missing from the assembly (black distribution between 10-40), which most likely represents haplotypes found only in the alternate phased assembly.

Figure 1

Distribution of kmers. One haplotype of the Brown Swiss genome compared to that of the complete read set (minus 10x barcodes). The red distribution signifies the distribution of read-kmers found only once in the assembly and the black distribution signifies the distribution of read-kmers not found in the assembly (mainly sequencing errors and haplotypes present in the alternate phased assembly).

We noted genes that were missing and duplicated in the Brown Swiss genome, ordering the duplicated genes in order of the most numerous duplications ( ). We note that all genes annotated as having five or more copies belonged to one of five genes. Three of these genes were repeat/transposon related genes. The remaining two genes were a predicted Bos indicus x Bos taurus elongation factor 1-alpha 1 pseudogene (6 copies), and the Bos taurus T cell receptor gamma cluster 1 (TCRG1) gene (128 copies).

Identification of Single Nucleotide Polymorphisms (SNPs) Within Bovine tlr2 Between Cattle Breeds

TLR2 is an important PRR involved in the recognition and uptake of MTBC family members. Key SNPs in TLR2 have been identified that affect susceptibility to mycobacterial infection for both humans and cattle (56, 57). The CDS of the bovine tlr2 gene of six BS and four HF cattle were cloned and compared to the bovine NCBI reference sequence genome (RefSeq) comprised of that of a B. taurus Hereford breed bull (NM_174197.2, NCBI, USA, accessed April 2017) in order to identify potential functional SNP sites. A total of 19 SNP variants were detected across the breeds, comprising both synonymous and missense mutations ( ). Of the eight missense mutations, seven occur in the ECD region of TLR2 and one in the TIR domain. One synonymous SNP at mRNA position 202 could not be classified with any available database, was described by Jann et al. (47) and is potentially novel. This position was called as ‘A’ in Water Buffalo, ‘G’ in BS and B. taurus, and ‘R’ (purine = G/A) in B. indicus, thus potentially indicating the presence of a heterozygous allele in B. indicus. Of all SNPs identified across the breeds, two of the non-synonymous SNPs have been already reported in the literature (44, 47, 56) and were only observed in the BS animals. One identified SNP variant, rs55617172, leading to an amino acid change of aspartic acid (D) to glutamic acid (E), both negatively charged, at amino acid position 63, has been significantly associated with TB resistance in a TB case- control cattle study (44). However, SNP rs68343167, was only identified in BS isolates in this study, and results in an amino acid change of the positively charged histidine (H) to the uncharged glutamine (Q) at amino acid position 326 (H326Q). This amino acid position has been suggested to be of relevance for ligand binding and functionality of human TLR2 (47, 58), bovine TLR2 and was additionally identified in the Anatolian Black Bos taurus breed (56). In addition, position 326 resides in leucine rich repeat 11 (LRR11) which together with LLR12 are the key domains involved in ligand biding and heterodimerisation of TLR2 with TLR1 and TLR6 (58).

Table 3

Identified SNP variants within the coding sequence of bovine tlr2 in HF and BS breeds animals.

								Brown Swiss						Holstein Friesian
dbSNP ID	Chromosome Position	mRNA Position	CDS Position	Amino Acid Position	Codon Change	Residue Change	Function of SNP	BS1	BS2	BS3	BS4	BS5	BS6	HF1	HF2	HF3	HF4
nd	nd	202	6	2	CCA->CCG	–	Synonymous	x	x	x	x	x	x	x	x	x	x
rs68268249	3953001	382	186	62	AAT->AAC	–	Synonymous		x	x	x	x	x
rs55617172	3952998	385	189	63	GAT->GAG	D(Asp) -> E(Glu)	Missense	x	x	x	x	x	x	x	x	x	x
rs68268250	3952985	398	202	68	GGC->AGC	G(Gly) -> S(Ser)	Missense		x	x	x	x	x
rs43706433	3952556	827	631	211	ATT->GTT	I(Ile) -> V(Val)	Missense	x	x	x	x	x	x	x	x	x	x
rs68268253	3970326	1141	945	315	CGG->CGT	–	Synonymous		x	x	x	x	x
rs68343167	3952209	1174	978	326	CAT->CAA	H(His) -> Q(Gln)	Missense		x	x	x	x	x
rs68343168	3952177	1206	1010	337	AGA->AAA	R(Arg) -> K(Lys)	Missense		x	x	x	x	x
rs68268256	3951937	1446	1250	417	AAC->AGC	N(Asn) -> S(Ser)	Missense		x	x	x	x
rs68268257	3951879	1504	1308	436	GGA->GGC	–	Synonymous		x	x	x	x
rs68268260	3951499	1884	1688	563	CGC->CAC	R(Arg) -> H(His)	Missense		x	x	x	x	x
rs440938496	3951477	1906	1708	570	GTG->GTT	–	Synonymous		x	x	x	x	x				x
rs68268261	3951408	1975	1779	593	GCA->GCC	–	Synonymous		x	x	x	x	x
rs467790258	3951405	1978	1782	594	GCG->GCA	–	Synonymous		x	x	x	x	x
rs68268263	3951192	2191	1995	665	CAC->CAG	H(His) -> Q(Gln)	Missense		x	x	x	x	x
rs68343171	3951162	2221	2025	675	CAT->CAC	–	Synonymous		x	x	x	x	x
rs68268264	3951132	2251	2055	685	ATT->ATC	–	Synonymous		x	x	x	x	x
rs68268266	3950973	2410	2214	738	GAG->GAA	–	Synonymous		x	x	x	x	x
rs68268267	3950892	2491	2295	765	CCC->CCT	–	Synonymous		x	x	x	x	x

The coding sequence for TLR2 in BS (n=6) and HF (n=4) animals were identified by Sanger sequencing and contigs matched against reference sequence NM_174197.2 (Bos taurus Hereford breed). SNP variants were classified as either synonymous or non-synonymous mutations, with the corresponding amino acid changes displayed. The SNP detected at mRNA position 202 was not found in the NCBI dbSNP database. The presence of allele change within individual animals is depicted as ‘x’. Amino Acid position highlighted depicts positions noted from TLR2 CDS transcript alignments detailed in . Yellow = unique to BS, Green = unique to BS and B. taurus, Blue = unique to B. indicus and Magenta = variable across all species.

Genotype Frequency Analysis of Candidate SNP H326Q in Selected Cattle Breeds

Having identified multiple SNPs within the coding sequence of tlr2 gene across breeds compared to the RefSeq genome, SNP rs68343167 (H326Q) was selected for further investigation using a larger sample size, including cDNA from two additional breeds, Boran and Sahiwal, from the B. indicus species. A total of 66 cDNA samples (Brown Swiss n=15, HF n=17, Sahiwal n=17, Boran n=17) were assessed for H326Q SNP analysis and frequency analysis ( ). The reference codon at the site of interest is encoded by the CAT sequence (TT genotype). Individual animals carrying the SNP had either heterozygote codons (TA genotype), encoded by either CAT or CAA, which is responsible for the H326Q variant; or they were homozygote for the CAA sequence (AA genotype).

Table 4

Genotype frequency of Bos taurus and Bos indicus cattle breeds for TLR2 selected candidate SNP rs68343167 (H326Q).

Cattle subspecies	Cattle breed	Sample size and (genotype frequency)			Total sample size (n)
		TT	TA	AA
Bos taurus	Holstein Friesian	17 (1.0)	0	0	17
	Brown Swiss	6 (0.4)	9 (0.6)	0	15
Bos indicus	Boran	2 (0.12)	8 (0.47)	7 (0.41)	17
	Sahiwal	0	0	17 (1.0)	17
	Total	25 (0.38)	17 (0.26)	24 (0.36)	66

A total of 66 cattle from Bos taurus and Bos indicus breeds (BS n=15, HF n=17, Sahiwal n=17, Boran n=17) were examined for their genotype frequencies with respect to SNP rs68343167 (H326Q). Genotype ‘TT’ codes for the homozygous ‘wildtype’ bovine TLR2 reference sequence NM_ 174197.2. Genotype ‘TA’ were heterozygous individuals and genotype ‘AA’ were homozygous for the H326Q SNP variant.

Interestingly, all HF samples tested as part of the present study were homozygous for the B. taurus Hereford RefSeq (TT genotype). Of the BS samples, six samples were homozygous for the RefSeq sequence, whereas the majority (n=9) were found to be heterozygous (TA genotype), containing both, the RefSeq sequence and the H326Q SNP. Of the Bos indicus breeds, the Boran samples revealed two samples homozygous for the RefSeq sequence, eight heterozygous and seven homozygous for the SNP. Thus, the heterozygote genotype was most abundant in the BS and Boran population. All Sahiwal samples (n=17) were homozygous for the SNP variant (AA genotype).

Alignment of TLR2 transcripts revealed additional SNP sites unique to BS, BS and B. taurus as well as B. indicus. Of note are positions H326Q as reported here and by others, as well as positions 417, 563 and 665; these latter 3 amino acid positions are particularly interesting as like H326Q these SNPs were found to introduce identical residue changes to our full length TLR2 analysis. Further, these SNPs were found to be present in BS cattle investigated ( ) but not present for HF cattle suggesting a role for genotype impact as found for H326Q. Position 417 and 665 are noteworthy, belonging to LRR15 and the TIR domain respectively. As LRR15 is close to the LRR regions thought to be involved in dimerization of TLR2 with co-receptors this SNP may play a role in effective PRR dimers forming. As the TIR domain is crucial for receptor dimerization of TLR cytoplasmic domains as well as signal potentiation; SNP H665Q may also modulate functional responses of TLR2.

On establishing distinct genotypes for the H326Q SNP in a range of cattle species (i.e., HF = TT, BS and Boran = TA and Sahiwal = AA), we assessed in the next step this association with TLR2 gene models and global sequence alignments. From the genotype analysis BS and Boran cattle represent an intermediate between the genotype extremes, perhaps representing genomic heritage; for example, Boran cattle are predominantly B. indicus with underlying influences of European and African B. taurus (59, 60). TLR2 gene models comparing the BS genome to B. taurus, B. indicus and B. bubalis (water buffalo) highlights an alignment of BS with B. taurus as TLR2 is derived from a single transcript in both breeds ( ). Highlighting the heterozygous genotype of BS cattle, BS resembled B. taurus at H326Q within LRR11 of the extracellular domain of TLR2 ( ) when the CDS of each possible TLR2 transcript was aligned to the three cattle species, Interestingly, additional sites where the BS and B. taurus sequences differed from B. indicus where I211V, F227L and R337K. Positions 211 and 337 were identified from the initial TLR2 sequencing described above.

The Brown Swiss H326Q TLR2 Polymorphism Induces a Stronger NF-κB Response After Ligand Binding Compared With Hereford and Holstein Friesian TLR2

Having identified the H326Q candidate SNP and its representation in the present study population of different cattle breeds, specifically in the BS population, the potential impact of this SNP on boTLR2 function in this breed was assessed. We utilised a previously published HEK cell assay to express the different TLR2 constructs in the absence of additional PRR (54). As well as the BS H326Q TLR2 variant, we included the Hereford reference variant and a L306P SNP variant that has been identified by others (47, 56) but not identified in the present study. We were also interested to assess the functionality of TLR2 across species, and therefore included the huTLR2 reference sequence (NCBI RefSeq Accession Number NG_016229.1). Furthermore, novel chimeric TLR2 constructs featuring the ECD of huTLR2, with the TIR domain of boTLR2 were created to delineate the contribution of each component to TLR2 functionality.

To confirm that the different TLR2 constructs were similarly expressed on the cell surface of the SEAP HEK cells, immunostaining was performed using bovine and human TLR2 specific antibodies. The mean fluorescence intensity (MFI) values of the positively stained population were similar for all constructs ( ). Stimulation with PMA, a non-specific activator of NF-κB, was used as a positive control and resulted in consistently high SEAP secretion for all transfectants. SEAP activity in response to ligands was normalised to PMA responses as an internal control.

The H326Q variant showed significantly higher TLR2-dependent NF-кB response in comparison to RefSeq Hereford TLR2 (p < 0.001 ). Constructs containing the L306P SNP did not respond significantly differently in comparison to the RefSeq Hereford boTLR2 version. The response of SEAP HEK cells transfected with huTLR2 was significantly higher in response to both FSL-1 and Pam3CSK4 in comparison to the bovine constructs. Interestingly, the chimeric bov:hu TLR2 receptor showed higher NF-κB activity than the RefSeq Hereford TLR2 HEK cell line, but significantly lower activity than the huTLR2 receptor in response to both FSL-1 and Pam3CSK4.

Figure 2

Location of identified SNPs. Gene model for TLR2 from assembly UMD 3.1.1 with annotated SNPs. Missense SNPs are located on the first track below the gene model, with synonymous SNPS on the second track. Reading from left to right, exons and coding sequence (CDS; thick line) are shown. All mutations occur in the second exon, which codes for ligand binding side.

NF-κB Activity Correlates Directly to CXCL8 Production of TLR2-HEK Cells Upon Stimulation With Lipopeptides

Using cell culture supernatant recovered from the above SEAP assay experiments, CXCL8 production in response to FSL-1 and Pam3CSK4 stimulation were measured ( ). As previously, there was an increase in CXCL8 production in the H326Q variant in comparison to the bovine RefSeq Hereford TLR2, however this was found only to be statistically significant in response to FSL-1 stimulation (p<0.01). CXCL8 production measured in huTLR2 transfectants was significantly higher than in boTLR2 transfectants in response to both ligands. CXCL8 production by chimeric bo:hu TLR2 HEK cells was significantly higher than any of the boTLR2 variants to FSL-1 (p<0.05), but lower than those from huTLR2, reflecting the observations made using the SEAP reporter assay of NF-κB activity.

Figure 3

TLR2-Ligand-Dependent NF-κB Activity and CXCL8 Secretion by HEK 293 cells Expressing Bovine and Human TLR2 Sequence Variants. HEK 293 cells, harboring NF-κB-induced SEAP reporter genes were transfected with either empty pTracer vector (mock); bovine or human TLR2 CDS constructs, or a bovine: human chimera in pDuo-mcs. Cells were stimulated with either FSL-1 at 100 ng mL-1 (A, C), Pam3CSK4 at 1 mg mL-1 (B, D) for 24 hrs. SEAP activity was measured by quantifying optical density of cell culture supernatant at 635 nm (A, B). OD values were normalised against values for PMA stimulation at 200 ng mL-1. Additionally, supernatants were assessed for CXCL8 concentration by ELISA (C, D). Error bars represent standard deviation from the mean of triplicate technical replicates and are representative of three independent repeats. (*=p<0.05, **=p<0.01, ***=p<0.001, ****=p<0.0001) and are only shown in relation to bovine TLR2 and additionally between the human and the chimeric bov: hu TLR 2 receptor. All data analysed by two-way ANOVA with multiple comparisons using GraphPad Prism V8 (GraphPad Inc., USA).

huTLR2 Is Significantly More Sensitive Than boTLR2 Variants to Stimulation With Live M. bovis BCG, but Not M. tuberculosis

Having assessed the activity of TLR2 variant HEK cells using TLR2-specific agonists and determining good correlation between SEAP activity and CXCL8 secretion ( ), we next challenged the TLR2-HEK cells with live mycobacteria and measured the NF-κB activation responses. Interestingly, the response of SEAP HEK cells expressing huTLR2 to challenge with M. tuberculosis was markedly lower than when they were challenged with TLR2-specific lipopeptides i.e., Pam3CSK4 and FSL-1 in the previous experiment ( ). Only the huTLR2 variant produced a significant NF-κB response upon challenge with M. bovis BCG compared to mock simulations. NF-κB responses for H326Q containing TLR2 were significantly lower than Hereford reference TLR2 and those constructs containing L306P. All boTLR2 variants responded more strongly to M. tuberculosis than to M. bovis BCG, even though the MOI was lower (MOI of 5) than with the M. bovis BCG infection (MOI of 10).

Brown Swiss and Holstein Friesian PBMC-Derived Macrophages Respond Significantly Differently to TLR2-Specific and Mycobacterial Ligands

Differences in functionality at the whole cell level conferred by boTLR2 breed variants were further assessed using PBMC-derived MØ. Cells were isolated from BS, containing the heterozygote H326Q variant, or HF cattle, which have the TLR2 sequence identical to the RefSeq Hereford breed. Due to Cat3 access restrictions, we measured the CXCL8 response to the TLR2-specfic agonist, FSL-1 and M.bovis BCG. In agreement with TLR2-HEK cell responses ( ), MØ generated from BS carrying the H326Q variant showed a significantly higher CXCL8 responses than HF MØ upon FSL-1 stimulation (p<0.05, ) and the 19 kDa protein, but not M.bovis BCG ( ). We further verified a high degree of correlation between the SEAP NF-kB and CXCL8 assays in a control experiment, by collecting data for both from the same stimulated cells and performing regression analysis ( ). Interestingly, unlike the SEAP HEK cell experiments, the PBMC-derived MØ responded similarly to MTB-derived 19 kDa lipoprotein and M. bovis BCG compared with FSL-1. We also detected a significantly greater response by BS derived MØ to M. bovis BCG compared to untreated cells, while the response by HF-derived MØ was not statistically greater than untreated cells.

Figure 4

M. bovis BCG and M. tuberculosis-Dependent NF-κB Activity Secretion by HEK 293 cells Expressing Bovine and Human TLR2 Sequence Variants. HEK 293 cells, harbouring NF-κB-induced SEAP reporter genes were transfected with either empty pTracer vector (mock); bovine or human TLR2 CDS constructs, or a bovine: human chimera in pDUO-mcs. Cells were stimulated with either live M. bovis BCG (MOI of 10) (A) or M. tuberculosis H37Rv (MOI of 5) (B) for 24 h. SEAP activity was measured by quantifying optical density of cell culture supernatant at 635 nm. OD values were normalised against values for PMA stimulation at 200 ng/ml. Error bars represent standard deviation from the mean of triplicate technical replicates and are representative of three independent repeats. (**=p<0.01, ***=p<0.001, ****=p<0.0001) and are only shown in relation to bovine TLR2 and additionally between the human and the chimeric bov: hu TLR 2 receptor. All data analysed by two-way ANOVA with multiple comparisons using GraphPad Prism V8 (GraphPad Inc., USA).

Figure 5

CXCL8 production by bovine MØ in response to FSL-1, M. bovis BCG or MTB Ligand Stimulation. Bovine MØ of BS (heterozygous for TLR2 H326Q SNP) and HF MØ (homozygous for the wild-type TLR2 sequence) were stimulated with either FSL-1 at 100 ng mL-1 (A), 19 kDa lipoprotein antigen (represents Rv3763 or Mb3789) (EMC microcollections, Germany) or infected with M. bovis BCG at MOI = 5 (B) for 24 hr. CXCL8 concentration in cell culture supernatant was assessed by ELISA. A total of n=8 animals per breed were stimulated with FSL-1 and n=6 animals per breed were infected with M. bovis BCG. Error bars represent standard deviation from the mean. Data was analysed with two-way ANOVA and Tukey’s post hoc comparison in GraphPad Prism V8 (GraphPad Inc., USA) (*=p<0.05, ***=p<0.001).

Discussion

Several studies have investigated cattle breed resistance to various diseases at the genetic level, leading to the suggestion that the HF breed may have lost immunological fitness, possibly due to extensive breeding for production traits (57, 61). Recently, we have also reported breed-specific differences in response to mycobacterial challenge (61). It is estimated that the cost of bovine tuberculosis to UK taxpayer is in excess of £100 million per year in surveillance testing and compensation, and there is currently no viable vaccine for cattle. Investigation into aspects of host immunity as they relate to different cattle breeds may inform future breeding strategies and minimise the impact of disease and reduce transmission of zoonotic TB to humans. Zoonotic TB is currently estimated to account for ~10% of global TB cases therefore taking a One Health approach to tackle zoonotic TB is paramount in achieving the UN Sustainable Development Goal of ending the TB epidemic (62). Reducing the burden of M. bovis infection in animal reservoirs such as livestock and wildlife underpins the 10 priority areas laid out in the WHO Roadmap for Zoonotic TB (62).

Using the power of linked-read sequencing, we have generated the first de novo phased assembly of the Brown Swiss cow genome. The genome shows high contiguity and completeness, as shown by the relatively high scaffold N50 and high single-copy ortholog reconstruction. Genome annotation lift over and gene duplication analyses show an increase in gene copy number of T cell receptor gamma cluster 1 (TCRG1) genes, in line with cattle belonging to “high γδ species” compared to human counterparts (63–65). It is important to note that γδ T cells are not only considered a major group of cells in mucosal immune responses, but also have been shown to be important in early immune responses to mycobacterial infection and bridge the gap between innate and adaptive immunity (66).

We identified 19 SNPs occurring within the CDS of bovine TLR2 compared to the RefSeq. A H326Q mutation was of particular interest as it resulted in an amino acid change from the positively charged histidine (H) to the uncharged glutamine (Q) at position 326, which is located in the ligand binding region LRR11 of the ECD of the receptor (58). Using a phylogeny-based approach, this SNP has been previously identified as a mutation under positive selective pressure in cattle and suggested to be of functional relevance (47). Furthermore, the H326Q variant was found by Bilgen et al. (56) in the B. taurus Anatolian Black cattle breed, which is thought to be more resistant to pathogens recognized by TLR2, such as to M. bovis and mastitis-causing bacteria. Among the B. indicus breeds, all Sahiwal samples in our study were homozygous for the H326Q SNP and except for two samples, half of the Boran samples were either homozygous or heterozygous for this variant, confirming a strong presence of this SNP in the B. indicus samples tested. This is interesting given that these breeds are known to be less extensively bred for production traits (56) and more resistant to several diseases, including bovine tuberculosis (57) compared with B. taurus. Within the BS breed, more than half were heterozygous for the SNP. For HF samples, the cattle breed under strongest selective pressure for production traits in our study population, all samples were homozygote and identical to the reference Hereford sequence. Since B. indicus cattle are phylogenetically more ancient than B. taurus cattle with respect to their common ancestor B. primigenius (67), our results suggest that BS cattle are immunologically potentially more closely related to ancient cattle breeds such as Boran and Sahiwal than the HF breed. Whether there is a co-evolution of different M. bovis spoligotypes with the bovine host is not yet know, but recent data suggest that it is more likely based on geographical area and farming practises rather than the host (68–70).

To compare the functional relevance of the TLR2 variants present in the HF and BS breeds, the identified sequences were cloned and expressed in HEK cells for in vitro phenotypic assessment. Initially, we used the SEAP assay that has been previously reported as a quantitative indicator of NF-κB activation during TLR stimulation (54, 71) and in which the SEAP response was directly correlated to CXCL8 production (22, 62). Upon stimulation with synthetic lipoproteins, constructs containing the H326Q SNP, showed a significant increase in NF-κB activity and therefore TLR2 signalling strength. As the H326Q SNP was well represented in our B. indicus population and most BS samples were heterozygous for this SNP, it can be hypothesized that this primes for a stronger TLR2-dependent immune response among these cattle.

For comparison, we included the characterisation of another reported SNP in the ECD domain of boTLR2, the L306P variant (47, 58) – although this was not identified in the present study. While this site is reported to be important for determining the size of the ligand-binding groove (58), we did not detect significantly different NF-kB or CXCL8 response to stimulation with TLR2-specific ligands used in this study. A downstream functional outcome of NF-kB activation is the transcription of certain inflammatory cytokines and including CXCL8 (72, 73), which is responsible for the recruitment of neutrophils to the site of infection and thus an important regulator of innate immunity. We observed a high degree of correlation between NF-kB activity and CXCL8 production in our assay, confirming the functional relevance of the SNPs between cattle breeds.

TLR function has been extensively studied, and shown to be species-specific and dependent on the ECD of the receptor (74). In contrast, there are conflicting reports over the contribution of the TIR domain, which is more highly conserved between species (7, 9, 74–76). We sought to clarify this aspect in the context of TLR2 expressed in isolation and challenged with TLR2-specific ligands. For both FSL-1 and Pam3CSK4, the bovine TIR showed effective, but reduced activity when expressed as a chimera with the huTLR2 ECD compared with the native human ECD-TIR structure. Similarly, the human ECD imparted greater TLR2 activity as a chimera with bovine TIR compared with native bovine ECD-TIR structure. Therefore, we can confirm that both ECD and TIR are independently responsible for differences in TLR2 sensitivity between species when TLR2 is expressed as a homodimer. In our study the human variant of both is more sensitive than the bovine orthologue to TLR2 ligands in the absence of other PRR. This contrasts with findings previously reported where the ECD was responsible for species-specific responses when a co-receptor construct of TLR2-Dectin1 was used (74), however these interactions may highlight the importance of co-receptor signal potentiation to increase the repertoire of innate immune responses. Our findings support the work of Faber et al. who demonstrated a similar result with human and porcine TLR5 chimeras (75).

Interestingly, the magnitude of the TLR2-dependent response was reduced when the cells were stimulated with live mycobacteria than with the synthetic ligands, independent of SNPs. There were some statistical differences between SNPs in NF-kB activity, but none of these were significantly different from the mock untreated conditions and so are not likely to be biologically relevant. Mycobacterium spp. have been reported to be capable of inhibiting the TLR2-driven response - for example, by the expression of glycolipids such as PDIM in the cell envelope that can obscure TLR2-ligand recognition (77). Furthermore, some MTB PAMPs may optimally require TLR2 heterodimerization or engagement of additional PRRs such as DC-SIGN to interact and induce NF-κB (14, 78, 79). By contrast, FSL-1 and Pam3CSK4 (80) have been shown to activate TLR2 as a homodimer (54) or TLR 2/6 heterodimers (81).

The response of SEAP HEK cells expressing bovine TLR2 was significantly higher when challenged with M. tuberculosis than with M. bovis BCG. Conversely, cells expressing huTLR2 were significantly more sensitive to M. bovis BCG than cells expressing boTLR2. As M. bovis BCG is derived from M. bovis, this may reflect species-level differences in host-pathogen coevolution that has resulted in humans being generally more resistant to M. bovis infection than cattle and vice-versa, cattle being more resistant to M. tuberculosis (82). M. bovis BCG is the avirulent vaccine strain of M. bovis that in certain circumstances induces protective immune but does not induce disease suggesting that modulated responses compared to pathogenic mycobacteria could be expected. It is interesting that this effect was observed when TLR2 was expressed in the absence of additional PRR. While heterodimerization between PRR has been shown to broaden the ligand spectrum for TLR2, it is not thought to induce differential downstream signalling pathways (83). The presence of accessory molecules such as CD14 and CD36 have been described to enhance TLR2 pathogen recognition by concentrating microbial products on the cell surface (84) and potentially this might be required for efficient boTLR2 recognition of M. bovis but is dispensable for recognition of M. tuberculosis. Further work is required to delineate the response of both species to virulent M. bovis in comparison with M. bovis BCG.

Whereas the HEK cell model enabled assessment of TLR2 function in isolation, we further characterised PBMC-derived MØ from either breed to examine a potential role in the context of additional PRR. BS MØ produced significantly higher CXCL8 to the TLR2-specific lipopeptide, FSL-1, than MØ from HF animals. Furthermore, BS MØ responded significantly better than unstimulated controls in response to M. bovis BCG, whereas HF MØ did not, and these data are in agreement with earlier work published by our group showing that M.bovis BCG does not induce a substantial response in MØ from HF animals (85). Taken together, our data are consistent with a positive role for the H326Q TLR2 variant supporting innate immunity at the whole cell level to mycobacterial infection. Another finding was that PBMC-derived MØ responded with greater CXCL8 production than HEK-TLR2 cells. This supports previously published observations that while TLR2 can initiate CXCL8 production, additional pathways can also upregulate the chemokine (11, 86). Therefore, the SNP in TLR2 contributes to, but is not solely responsible for the inflammatory response to mycobacteria. It is likely that additional breed-specific factors, including a possible contribution by adaptive immunity determines the bovine response to tuberculosis at the whole animal level.

Our findings provide proof of concept however, that meaningful genetic differences exist between cattle breeds that impact immune function, and such work can be expanded to consider additional factors. It is also interesting to speculate which of the identified TLR2 sequences may actually represent the “wild-type” sequence and whether what we see in cells from HF is due to their intensive genetic selection breeding for milk-production, resulting in an immunologically less-fitter status compared to the sequence/responses seen in cells from BS animals. Resistance to disease is often multifactorial and polygenic and thus identified resistance markers are often unlikely to provide absolute resistance to disease (87). Of course selection for resistance should not occur at the expense of other control measures such as biosecurity and bTB control in the wildlife reservoir (88). Immune function and mycobacterial resistance factors associated between cattle breeds (and species) should be considered when implementing vaccine-based or immunomodulation-based control measures to reduce M. bovis burden. However, more importantly, the data presented in our work, as well as in the work published by authors already point towards the possibility that in the long-term, breeding for bTB resistance could definitively be achieved. We do not envisage that TLR2 will be the only PRR under selection, but indeed may be one PRR worth further investigation. Whether this increase in disease resistance would go hand-in-hand with a potential beneficial reduction in production, but increase in animal longevity and thus fewer replacements remains to be seen,

Data Availability Statement

The datasets presented in this study can be found in online repositories. The Brown Swiss genome assembly, along with the associated raw read data are available from the European Nucleotide Archive under study accession PRJEB47135.

Ethics Statement

All procedures were carried out under the Home Office license (PPL7009059) approved by the RVC’s AWERB Committee.

Author Contributions

Writing proposal: AG, SW, and DW. Experimental work and scientific conceptualization: M-CB, AG, SW, AH, and DW. Genome assembly, analyses, and annotation: GE and FP. All authors contributed to the article and approved the submitted version.

Funding

The project was funded by a Bloomsbury PhD studentship to AG, SW and DW, as well as BBSRC grant BB/P008461/1 to DW. The manuscript is ID No 1550801 of the RVC.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.