Genome-wide association analysis of avian resistance to Campylobacter jejuni colonization identifies risk locus spanning the CDH13 gene

Abstract

The enteropathogen Campylobacter jejuni is a major worldwide health and economic burden, being one of the leading causes of bacterial gastroenteritis and commonly linked to postinfectious onset of autoimmune disease. Chickens are a major vector for human infection and even though variation in avian colonization level is heritable, no previous studies have identified regions of the genome associated with colonization resistance. We performed a genome-wide association study of resistance to C. jejuni colonization in the avian intestine by controlling for population structure, which revealed a risk locus with genome-wide significance spanning the T-cadherin (CDH13) gene. A second possible risk locus was also identified close to calmodulin (CALM1), a calcium-activated modulator of cadherin function. In addition, gene expression analysis of mRNA sequencing profiles revealed that the relative expression of the two genes is significantly associated with colonization resistance. Functional studies have previously demonstrated involvement of cadherins and calmodulin in C. jejuni intracellular invasion and colonization of human intestinal epithelial cells in vitro. Consistent with this finding, our analysis reveals that variation surrounding these genes is associated with avian colonization resistance in vivo and highlights their potential as possible targets for control of the bacterium in avian and human populations.

Keywords: GWAS; cadherin; calmodulin; epithelial cell invasion; immunity; intestinal homeostasis.

Figure 1

Colonization status of chickens selected for the genome-wide association study. (A) The initial population of 255 chickens had highly variable colonization levels; 196 of these were selected for genotyping using the 60K chicken SNP chip. Their colonization levels are illustrated in (B). A total of 168 chickens from opposite ends of the colonization spectrum (i.e., all individuals with colonization level ≤400 CFU/g and all individuals with colonization level ≥10⁷ CFU/g) were selected to conduct the case-control association study. Two individuals were filtered due to genotyping call-rate <95%. Colonization levels of all 166 individuals which passed quality control filtering for the case-control study are plotted in (C). Forty of these had no or very minimal colonization (resistant cases), and 126 had high colonization levels (susceptible controls). Birds for which RNAseq expression profiles are available are highlighted in red at the extremes of the distribution. The lower limit of detection = 400 CFU/g.

Figure 2

Principal components revealing population stratification. All 194 individuals that passed quality control filtering are represented. Colonization-susceptible individuals from the case-control study are shaded blue and colonization-resistant birds shaded red. Gray circles represent individuals from the center of the colonization spectrum that were not used in the association study. All autosomal markers were used to construct a kinship matrix of genetic relatedness or identity by descent (IBD). This matrix was converted to a distance matrix and classical multidimensional scaling performed to detect any population substructure. There is evidence for some stratification of the population which was resolved by the first three principal components (PC1, PC2, and PC3). The percentage variation which is represented by each of these principal components is 18.6% (PC1), 16.8% (PC2), and 11.1% (PC3). The first two principal components reveal separation of the population into two subpopulations—one large group with 157 individuals and one smaller group of 37 individuals.

Figure 3

Quantile-Quantile plot of observed vs. expected P-values. The q-q plot of observed P-values shows they follow the expected null distribution apart from 12 SNPs with the most significant P-values. These are located within three regions: chromosomes 11 (red squares), 5 (blue triangles), and 1 (yellow circles).

Figure 4

Manhattan plot of genome-wide association. All 16,871 SNPs that passed quality control are shown. Genomic location is plotted against −log₁₀(P). Four SNPs on chromosome 11 reached genome-wide significance (FDR = 0.04).

Figure 5

Regional plot of significant association on chromosome 11. Four SNPs with genome-wide significance FDR < 0.05 are highlighted red. These four SNPs are intergenic but surround the T-cadherin (CDH13) gene. Two of these SNPs are in complete LD, situated very close to each other, and are almost indistinguishable on the plot. The physical positions of all Refseq genes within the region are depicted below the plot.

Figure 6

Relative expression levels of CDH13 and CALM1. Boxplots of resistant vs. susceptible birds reveal separation of the population based on CDH13/CALM1 expression (Mann-Whitney U-test, P = 0.004) and no bias in relative ratios due to sequencing depth. The ratio of CDH13/CALM1 expression for all 28 transcriptionally profiled birds is plotted, with shading reflecting depth of sequencing.