Multiple apical plasma membrane constituents are associated with susceptibility to meconium ileus in individuals with cystic fibrosis

Variants associated with meconium ileus in cystic fibrosis (CF) were identified in 3,763 patients by GWAS. Five SNPs at two loci near SLC6A14 (min P=1.28×10−12 at rs3788766), chr Xq23-24 and SLC26A9 (min P=9.88×10−9 at rs4077468), chr 1q32.1 accounted for ~5% of the phenotypic variability, and were replicated in an independent patient collection (n=2,372; P=0.001 and 0.0001 respectively). By incorporating that disease-causing mutations in CFTR alter electrolyte and fluid flux across epithelia into an hypothesis-driven genome-wide analysis (GWAS-HD), we identified the same SLC6A14 and SLC26A9 associated SNPs, while establishing evidence for the involvement of SNPs in a third solute carrier gene, SLC9A3. In addition, GWAS-HD provided evidence of association between meconium ileus and multiple constituents of the apical plasma membrane where CFTR resides (P=0.0002, testing 155 apical genes jointly and replicated, P=0.022). These findings suggest that modulating activities of apical membrane constituents could complement current therapeutic paradigms for cystic fibrosis.

P=0.022). These findings suggest that modulating activities of apical membrane constituents could complement current therapeutic paradigms for cystic fibrosis.
Cystic fibrosis is caused by mutations in the cystic fibrosis transmembrane conductance regulator (CFTR) gene 1 . CFTR is a chloride channel located on the apical membrane of epithelia, where ion conduction and solute trafficking contribute to the regulation of transepithelial fluid flow. Individuals with the same loss-of-function CFTR mutations have variable disease severity, and differentially affected CF-associated organs including lung, pancreas, liver, intestines, and vas deferens; thus additional features, including other genes (referred to as modifier genes) may affect disease pathophysiology. Approximately 15% of CF patients have severe intestinal obstruction at birth, a complication known as meconium ileus 2 . Meconium ileus develops in utero, and presents following birth with complete intestinal obstruction that requires either medical or surgical intervention. This neonatal complication is highly indicative of CF, occurs in either sex, displays notable heritability exceeding 88% 3 , and is likely minimally affected by environmental influences.
The North American CF Gene Modifier Consortium has accumulated 3,763 participants with 'severe' (pancreatic exocrine insufficient) CFTR alleles and genome-wide genotype data at 543,927 SNP loci 4,5 (Table 1 and Online Methods). The definition of meconium ileus was consistent within the consortium and was recorded following rigorous chart review. A conventional GWAS for meconium ileus used a generalized estimating equations (GEE) model 6 to include collected sibling-pairs, and led to five genome-wide significant SNPs (P<5×10 −8 ) 7 from two regions that include SLC26A9 on chromosome 1 and SLC6A14 on chromosome X (Fig. 1, Supplementary Fig. 1, Table 2; sex-specific results in Supplementary Table 1). CFTR was not a significant confounder or effect modifier when incorporated in the GWAS (Supplementary Fig. 2 and Supplementary Table 2), indicating SLC6A14 and SLC26A9 are independent contributors to meconium ileus. We then replicated the associations in SLC6A14 (min P=0.001) and SLC26A9 (min P=0.0001) with meconium ileus in an independent combined collection from North America and France ( Table 2).
The signal intensity plots of the associated SNPs reflected autosomal-and X-associated SNPs at SLC26A9 and SLC6A14, respectively. Imputation analysis using MACH and minimac 8,9 identified the same regions of association as the genotyped SNPs (Online Methods, Supplementary Fig. 3). The five associated SNPs in SLC6A14 and two in SLC26A9 ( Fig. 1b and 1c) are positioned just upstream of their respective transcription start sites such that binding of activating or repressing transcription factors may be affected as highlighted by ENCODE data 10 (data not shown). Neither SLC6A14 nor SLC26A9 coding regions exhibit evidence for CNVs; however, there is a gap in the reference sequence >10 kb upstream of the SLC26A9 locus.
The seven SNPs genotyped (Table 2) in the two genes account for <5% of the meconium ileus variation, estimated by pseudo R-squared 11 , likely reflecting the common problems in association studies of locus heterogeneity and low power given the available sample size. Whereas conventional GWAS is often considered for complex disease mapping, modifier gene studies could incorporate disease etiology and pathobiology information to increase power and account for heterogeneity. To do so, we proposed application of GWAS with consideration of a hypothesis (and so is hypothesis-driven; GWAS-HD) to systematically prioritize SNPs for genome-wide analysis. The highest priority markers are also evaluated as a set to test the statistical significance of the hypothesis used for prioritization ( Supplementary Fig. 4).
The GWAS-HD prioritization in this CF application is based on the knowledge that a major source of CF pathophysiology is impaired fluid and electrolyte flux in epithelia of CFaffected organs. The polarized epithelial layer forms a highly selective barrier between organ and ductal interfaces. Transepithelial 'function' is achieved by cell polarization whereby many determinants and regulators of fluid, solute and ion transport reside at the apical membrane alongside CFTR, with contributing features from basolateral surfaces. We have shown in a mouse model that CFTR function in the gastrointestinal epithelium is critical for preventing intestinal obstructions 12 . Thus, we hypothesized that with loss of CFTR, (genetic) variation in other apical membrane constituents could modify CF phenotypes, such as meconium ileus.
A list of 157 gene products ( Fig. 2 Table 3) was annotated as localized to the apical plasma membrane using AmiGO 13 with Gene Ontology data 14,15 . CFTR and many solute transporters were included. However, the brush border membrane protein SLC6A14 was not listed, reflecting the high specialization of its corresponding intestinal cavity and a limitation of the GO annotation that we accepted without additional curation to avoid bias. In total, 3,814 GWAS SNPs were within ±10 kb of the boundaries of 155 genes (NCBI36/hg18); 2 genes were not tagged by any of the GWAS SNPs.

and Supplementary
To implement the GWAS-HD for meconium ileus using the apical hypothesis, we first prioritized the genome-wide markers by assigning the 3,814 SNPs of the apical genes to a high priority group and all remaining genome-wide SNPs to a low priority group. We then performed two statistical analyses (Supplementary Fig. 4 and Online Methods). The first, analogous to a conventional GWAS, was to conduct single-SNP association analysis using all of the 543,927 GWAS SNPs at the genome-wide level, however after up-weighting the 3,814 apical SNPs via the SFDR control. 16 The second analysis, focusing only on the 3,814 high priority SNPs using a multi-SNP/gene analysis, tested the prioritization hypothesis itself to determine if multiple proteins present on the apical plasma membrane contribute to meconium ileus susceptibility.
As in the conventional GWAS, SNPs from SLC6A14 showed the strongest evidence for association with meconium ileus in the single-SNP GWAS-HD analysis ( Supplementary  Fig. 5), despite not being in the high priority group, reflecting the robustness of SFDR 17 . In addition, SNPs from SLC26A9, and two additional apical genes, ATP2B2 and SLC9A3, showed association evidence with q values <0.05 (Table 3). A gene-based analysis (Online Methods) of ATP2B2 (P=0.0006) and SLC9A3 (P=0.0001) indicated evidence for allelic heterogeneity after comparing results with single-SNP analysis (Table 3). SLC9A3 was replicated in the French cohort (P=0.017) while ATP2B2 was not (P=0.283) (Supplementary Table 3).
Next, restricting analysis to the 3,814 SNPs annotated to the 155 apical genes (which does not include SLC6A14), we tested the apical prioritization hypothesis as part of the GWAS-HD. Here we observed genome-wide significant evidence for association between meconium ileus and multiple constituents of the apical plasma membrane (permutation P=0.0002, testing all 3,814 SNPs jointly and not subject to multiple hypothesis testing ( Fig.  2a and 2b)). Even with the exclusion of SLC26A9 (as well as SLC6A14), the apical hypothesis remained significant (P=0.0058). Thus, GWAS-HD further established the involvement of other genes coding for apical constituents despite insufficient power to detect individual SNPs, even within the context of our prioritized GWAS. For comparison, we also constructed a null hypothesis list of membrane-localized genes. As expected, the 224 GO-annotated nuclear envelope genes tagged by 3,537 GWAS SNPs showed no relationship with meconium ileus (permutation P=0.4639; Fig. 3).
The French cohort with genome-wide data provided independent validation of the apical hypothesis (permutation P=0.022; Fig. 2c and 2d; Online Methods). The statistical significance of this gene set (which excludes SLC6A14) in the French cohort remained after further excluding SLC26A9 (P=0.021) and then both SLC26A9 and SLC9A3 (P=0.023). Although analysis of the apical hypothesis in a larger independent cohort should be considered as part of future efforts, the replication in the smaller French cohort supports a (common) mechanism of contributing genes, even though the responsible gene variants across the two datasets may not be the same.
To determine which apical genes were driving the association, the degree of genetic heterogeneity in meconium ileus, and the common contributors across the French and North American samples, we implemented Lasso 18 . Using the North American sample, we jointly analyzed all 3,740 SNPs available in the apical genes (which include SLC26A9 and SLC9A3), and SLC6A14 (Online Methods), that were not in perfect linkage disequilibrium. Forty-eight SNPs spanning 36 different genes were retained by Lasso in the multivariate regression model (Supplementary Table 3). These SNPs jointly explained ~17% of the meconium ileus variation in the North American sample. The percentage explained by the same 48 SNPs in the French sample was 8.1% (Online Methods), presumably due to the smaller sample size and genetic heterogeneity. We then tested the significance of a score for each individual in the French cohort constructed from a weighted sum of the number of risk alleles (defined in the North American sample) across the 48 SNPs 19 . The significant association between meconium ileus and this score (P=0.0137; Online Methods) provided additional complimentary evidence of common contributors between the two cohorts, with SNPs in SLC9A3 and SLC6A14 being two specific examples.
In summary, conventional GWAS identified SNPs in SLC6A14 and SLC26A9 as significantly associated with meconium ileus. GWAS-HD single-SNP analysis identified the same SNPs in SLC6A14 and SLC26A9, as well as SNPs in SLC9A3; and multi-SNP analysis provided evidence that multiple constituents of the apical plasma membrane are collectively associated with meconium ileus, yielding considerable additional information beyond single SNP/gene associations. GWAS-HD can be applied to other Mendelian disorders, or even complex traits provided there is a biologically-based hypothesis and participating relevant genes can be compiled.
Although gene prioritization has been used in other approaches such as pathway or gene enrichment analyses [20][21][22] , GWAS-HD involves key differences. First, in contrast to the previous inclusion/exclusion approaches where all genotyped SNPs/genes are not analyzed simultaneously, GWAS-HD performs parallel single-SNP analysis of all GWAS SNPs, and multi-SNP/gene analysis focusing on the set of SNPs/genes of interest. The prioritized single-SNP analysis interrogates all available SNPs via the SFDR control, yet enables increased statistical power for regions favoured a priori. For example, SLC6A14 would be omitted by inclusion/exclusion approaches, yet it remained the highest ranked gene for association with meconium ileus in GWAS-HD. Second, methods such as interactive pathway analysis can be restrictive because contributing genes/proteins must relate to each other via direct or indirect links, which may be disturbed when a cog component (such as CFTR) is dysfunctional in the disease state. While gene products that participate in maturation or delivery of CFTR may be contributory, consideration of local components of processes that may compensate for the ion and fluid flow disturbance in CF is enabled in the apical hypothesis. Third, distinct from an exhaustive search of all plausible interactive pathways, GWAS-HD focuses on a single biological hypothesis and provides statistical significance for all genes involved jointly, alleviating the multiple testing burden.
It should be noted that the specific statistics or models used in our GWAS-HD application such as SFDR, Lasso, and the sum and score statistics may not be the most powerful ones in any specific setting. For example, the adaptive rank truncated product statistic 23 could be used to identify the common subset of associated apical genes across two samples; and there are alternative weighting and prioritization approaches 24 .
The SLC9A3 gene codes for a sodium/hydrogen exchanger that when disrupted has been shown to decrease intestinal obstructions in a CF mouse model 25 . SLC6A14 codes for a sodium-and chloride-dependant neutral and basic amino acid transporter 26,27 . SLC26A9 encodes an anion transporter, likely a chloride channel with multiple modes that include chloride/bicarbonate exchanger and sodium-anion cotransporter capabilities 28 . It has also been reported to physically interact with CFTR 29 and be influenced by CFTR activity, at least in lung-related tissues 30 .
It is notable that SLC9A3 has previously been associated with infections and pulmonary function in CF 31 . In the consortium discovery sample, rs6864158 (MAF=0.43) in SLC9A3 was associated with both the lung (P=0.0003, analyzed previously 5 ) and meconium ileus (P=0.0001); this provides evidence that some genes may play a role in multiple CF comorbidities. Both SLC6A14 and SLC26A9 are robustly expressed in human lung epithelia, sweat gland, as well as intestinal epithelia as measured by RT-PCR (not shown). We anticipate that meconium ileus modifier genes may also influence early pathology in other CF-affected organs, and thus could provide important insights into the mechanisms of CF disease severity and co-morbidities.  Meconium Ileus GWAS identifies genome-wide significant SNPs. Association analysis was performed on all SNPs with minor allele frequencies > 2% that passed QC criteria (Online Methods). (a) Genome-wide Manhattan plot for meconium ileus. The black solid line corresponds to the genome-wide significance threshold 7 (P<5×10 −8 ), and the black dashed line to the suggestive association threshold, expected once per genome scan (P <1/543,927=1.84×10 −6 ). A total of five SNPs in two regions (SLC6A14 on chromosome X, and SLC26A9 on chromosome 1) have association evidence exceeding the genome-wide threshold. The SNPs, rs4077468 and rs4077469 are in perfect LD and appear as one SNP as they are separated by only 128bp. (b) Regional plot for SLC26A9. LocusZoom viewer was used to display the association evidence around SLC26A9 based on NCBI Build 36/hg18. Symbol coloring reflects HapMap CEU LD r 2 values with the most significant SNP. The significant SNPs, rs4077468, rs7419153 and rs12047830 (Table 2), occur 2.17 kb, 4.72 kb and 4.12 kb upstream, respectively, of the transcription start site. The significant SNP, rs7512462, occurs in intron 5. A gap occurs in the genomic sequence between the SLC26A9 and FAM72A genes in both NCBI36.3 and GRCh37 primary reference assemblies.
(c) Regional plot for SLC6A14. The association evidence around SLC6A14 was displayed as above. Rs3788766 (Table 2) is positioned 0.95 kb upstream of the transcription start site and is within the binding site of the CEBPB transcription factor as annotated by ENCODE 33 (not shown). The mRNA transcript corresponding to CXorf61 is downstream (3') of SLC6A14. The apical membrane hypothesis identifies genes associated with meconium ileus. A list of 157 genes was annotated using the AmiGO tool 13 (version 1.7; March 28, 2010) based on the Gene Ontology data 14 (GO:00163245) using the cell location search phrase "apical plasma membrane" with restriction to Homo sapiens. In total, 3,814 GWAS SNPs are within ±10 kb of the boundaries of 155 genes (NCBI36/hg18). Two genes were not tagged by any of the GWAS SNPs; SLC6A14 was not annotated to the apical plasma membrane. (a) QQ-plot of the apical SNPs in the discovery sample. The observed association statistics (red), and the statistics calculated from the 10,000 phenotype-permutated replicates are shown (light gray). (b) Statistical significance of the apical membrane hypothesis in the discovery sample. Statistical significance (permutation P = 0.0002) was established via a sum statistic, summing the association evidence (Wald χ 2 1 statistic) over all the 3,814 SNPs with the Assessment of the nuclear envelope null hypothesis. A list of 231 genes was generated from the nuclear envelope as defined by GO annotation (GO:0005635) similarly as for the apical membrane list. In total, 3,537 GWAS SNPs are within ±10 kb of the boundaries of 224 tagged genes (NCBI36/hg18). A priori, the nuclear envelope list should not contain genes associated with meconium ileus (under the null of no association). (a) QQ-plot of the nuclear envelope gene SNPs in the discovery sample. (b) Statistical significance of the nuclear envelope hypothesis in the discovery sample. Statistical evaluation indicates that genes listed in the nuclear envelope are, as expected, not significantly associated with meconium ileus (permutation P=4639/10,000=0.4639). Sun Table 3 Ranked SNPs with q values < 0.05 from GWAS or GWAS-HD

SLC9A3 a
The q value is a genome-wide adjusted P value that controls the false discovery rate. b The rank indicates the genome-wide ordering of the SNPs based on the original association evidence (GWAS) or SFDR q value after incorporating the apical plasma membrane hypothesis (GWAS-HD).
The GWAS and GWAS-HD rank results of all 543,927 SNPs are provided in Supplementary Fig. 5. c Apical indicates whether a gene was (= 1) or was not (= 0) on the generated apical membrane constituent list.