Identification of three new cis-regulatory IRF5 polymorphisms: in vitro studies

Genotypes, linkage disequilibrium map and haplotype analysis

The region extending from chr7:128,343,695-128,393,203 of the Human March 2006 genome assembly (NCBI36/hg18) that includes the coding sequence of IRF5, 21.5 kb of its promoter region (from the start of exon 1, variant A) and part of the transportin 3 (TNPO3) locus, which is 3' to IRF5, was screened for polymorphisms. Only tag SNPs were selected in about half of this region; in the promoter and first intron of IRF5, however, we selected all SNPs characterized by double-hits in NCBI dbSNP Build 128 [18]. The 54 selected SNPs (Additional file 2: Table S1) were analyzed in DNA samples obtained from peripheral blood mononuclear cells (PBMCs) of 95 Spanish healthy controls recruited in our center. They provided their written informed consent to participate in the study, and the study was approved by the Ethics Committee for Clinical Research of Galicia. Genotypes were obtained by single-base extension with the SNaPshot Multiplex Kit (Applied Biosystems, Carlsbad, CA, USA) as described previously [5], except for rs3778752, rs3778751 and the CGGGG indel, which were sequenced; the exon 6 indel, which was genotyped by length variation in agarose gel electrophoresis as described previously [5]; and rs10954213, which was genotyped using a fluorogenic 5' nuclease assay (TaqMan MBG Probes, TaqMan SNP Genotyping Assay; Applied Biosystems). The primers and probes used are given in Additional file 2: Table S2. The genotyping call rate was 97.7%. All the SNPs were in Hardy-Weinberg equilibrium with a threshold for significance of 0.05 without correction for multiple testing. Pairwise D' values, r² values and their graphic representation were obtained using Haploview software at the default settings [19]. Haplotypes were estimated using PHASE 2.1 software (at default parameters), which implements a very accurate Bayesian algorithm [20]. Two haplotype distribution analyses were done, the first to select SNPs for functional studies, including all SNPs of potential interest and the second to put the new functional polymorphisms in context, including only the functional and tag SNPs. Functional screening was conducted on three additional polymorphisms selected on the basis of their relationship with IRF5 expression or their known cis-regulatory effect. This information is detailed in Additional file 2: Table S3.

Cell culture and nuclear extracts

The B-lymphocyte cell line WIL2 NS [21, 22], which constitutively expresses IRF5 (Additional file 3: Figure S1), was obtained from the European Collection of Cell Cultures (catalogue no. 90112121). It was maintained in RMPI 1640 medium supplemented with 10% fetal bovine serum, 2 mM L-glutamine, 1 U/ml penicillin and 0.1 mg/ml streptomycin. For the preparation of nuclear extracts, 10⁷ cells were washed in 1 ml of phosphate-buffered saline and incubated at 4°C for 15 min with 200 μl of a lysis buffer containing 10 mM 2-[4-(2-hydroxyethyl)piperazin-1-yl]ethanesulfonic acid (HEPES), pH 7.9, 1 mM ethylenediaminetetraacetic acid (EDTA), 1 mM ethylene glycol tetraacetic acid, 10 mM KCl, 1 mM dithiothreitol (DTT), 1 mM phenylmethylsulfonyl fluoride, 10 μg/ml aprotinin and 10 μg/ml of leupeptin. The freed nuclei were incubated at 4°C for 30 min with 0.1% Triton X-100 and collected by centrifugation at 800 × g during 15 min at 4°C. Subsequently, they were incubated at 4°C for 30 min with rotation in 200 μl of the previously specified lysis buffer supplemented with 20% glycerol and 0.4 M KCl. Insoluble material was precipitated by centrifugation for 15 min at 13,000 rpm and 4°C. All the mentioned reagents were obtained from Sigma-Aldrich (St Louis, MO, USA). Nuclear proteins in the supernatant were quantified using the Quant-iT Protein Assay Kit (Molecular Probes, Eugene, OR, USA) in a Qubit fluorometer (Invitrogen, Carlsbad, CA, USA). Quantified nuclear extracts were stored at -80°C until use.

Electrophoretic mobility shift assay

Oligonucleotides for each allele and strand centered in the 25 selected polymorphisms plus the CGGGG indel used as a positive control were synthesized by Sigma-Genosys (Haverhill, UK) (Additional file 2: Table S4). Complementary oligonucleotides were annealed at 37°C for 1 h to form the double-stranded probes used in EMSA. The LightShift Chemiluminescent EMSA Kit (Pierce Biotechnology, Rockford, IL, USA) was used with some modifications. Briefly, the binding reaction was performed with 0.5 nM of double-stranded, biotin-labeled oligonucleotides and 0.5 to 4.0 μg of nuclear extract proteins from WIL2 NS cells. This reaction was brought to a final volume of 20 μl with binding buffer (10 mM Tris·HCl, pH 7.5, 10 mM KCl, 0.6 mM DTT, 50 ng/μl poly(deoxyinosinic-deoxycytidylic) acid, 1.5 mM EDTA, 2 mM HEPES, pH 7.9, 25 mM NaCl, 0.1 mM ZnSO₄, 15% glycerol and 0.25 mg/ml bovine serum albumin) and incubated at room temperature for 20 min. Cold competition was done with a 200-fold molar excess (0.1 μM) of the unlabeled oligonucleotide. DNA-protein complexes were analyzed by electrophoresis in 5% polyacrylamide nondenaturing gels and 0.5 × Tris/borate/EDTA buffer at 11 V/cm and 4°C for 3 h. Subsequently, they were transferred to nylon membranes according to a semidry transfer protocol for 45 min at 12 V. The complexes were cross-linked to the nylon membrane with ultraviolet (UV) light at 7.2 J/cm². Detection was done with the LightShift Chemiluminescent EMSA Detection Module (Pierce Biotechnology) followed by exposure to UV light using the UVP EC3 Imaging System (UVP, Upland, CA, USA). A minimum of three independent replicate experiments were done for each EMSA showing differential binding.

Luciferase reporter assays

Fourteen different constructs were made, one for each of the alleles of the seven polymorphisms with differential EMSA results. The IRF5 sequence (148 to 245 bp in length) of each of these polymorphisms was obtained from heterozygous subjects by polymerase chain reaction (PCR) with the same primers used for genotyping (Additional file 2: Table S5). None of them included additional polymorphisms. The amplicons were first inserted into the EcoRV site of a pBlueScript (pBK SK-) vector by TA cloning. Briefly, the EcoRV cut vector was treated with Taq DNA polymerase in the presence of 2 mM deoxythymidine triphosphate and PCR buffer, and 50 μg of this 3' T-vector were ligated to 25 μg of the PCR amplicon in a 10-μl reaction with 1 μl of T4 DNA ligase reaction buffer (New England Biolabs, Ipswich, MA, USA) and 1× ligase buffer at 16°C overnight. The cloned IRF5 sequences were subcloned into the KpnI/SmaI sites of a fos-pGL3-basic vector, which are just upstream of the minimal fos promoter regulating expression of the Firefly luciferase gene. The sequence and orientation of each insert were verified by DNA sequencing. These verified reporter vectors were transfected in WIL2 NS cells. Specifically, each experiment included five conditions: blank, positive control, empty vector, allele 1 and allele 2. The experiments were carried out by transfecting 2 × 10⁶ WIL2 NS cells with 400 ng of the adenovirus E1 region RSV-Luc vector plus 2,800 ng of the pBK SK- vector for the positive control, 3,200 ng of the fos-pGL3 basic vector for the empty vector, 3,200 ng of the allele 1 fos-pGL3 vector construct, and 3,200 ng of the allele 2 fos-pGL3 vector construct for each of the two alleles. All transfections also included 800 ng of the pRL-TK Renilla luciferase vector (Promega, Madison, WI, USA) for normalization. Transfections were done by microporation in 100 μl of Resuspension Buffer R at 1,100 V for 30 ms and a single-pulse program in the Neon Transfection System (Invitrogen). After microporation, the transfected cells were left in culture for 24 h (basal) or left for 20 h. The luciferase activity was measured using the Firefly Luciferase and the Renilla Luciferase Assay systems (both from Promega). The results for Firefly luciferase were expressed as relative luciferase units after normalization with the Renilla luciferase activity and with the empty vector for the respective experiments. At least five independent experiments were performed for each of the seven analyzed polymorphisms. Comparison of the log-transformed mean normalized Firefly luciferase luminescence ratios was done using the Wilcoxon matched-pairs test.

Selection of promoter polymorphisms on the basis of linkage disequilibrium and haplotypes

The experimental design involved selection of SNPs on the basis of three criteria (Figure 1): high r² values with rs729302, high correlation with the protective haplotypes and prominent association with IRF5 expression. These three criteria are not mutually exclusive, and many SNPs were selected for more than one of them. The selected SNPs were evaluated for their effect in binding nuclear proteins as an indication of their possible role in cis-regulation of IRF5 expression by means of EMSA. The SNPs with differential binding to their alleles in EMSA were assessed with reporter gene assays, and those showing changes in expression between the constructs with the two alleles were considered to be functional (Figure 1).

To identify SNPs strongly correlated with rs729302 or with the protective haplotypes, we genotyped 54 polymorphisms located predominantly within 21.5 kb of the 5' region of IRF5 (Additional file 2: Table S1) in 95 Spanish controls. The D' and r² pairwise values were calculated for 44 polymorphisms (Additional file 2: Table S6), which excluded eight SNPs with minor allele frequencies less than 2% (rs12536195, rs6951615, rs11773414, rs11983607, rs11763323, rs754280, rs41298401 and rs11767834) and two SNPs that were completely redundant with rs3778753 (rs3778752 and rs3778751). This analysis showed a block of LD extending over 43 kb from TPNO3 to the 5' region of IRF5 further from the gene than rs729302 (Figure 2A). This block included 37 of the analyzed SNPs. A second block of LD with seven SNPs was even further 5' to IRF5 than the first. In the largest block of LD, there were two areas of stronger correlation between SNPs (Figure 3B), one of them (8 kb) included rs729302 and other eight SNPs (r² ≥ 0.7 with rs729302). These nine SNPs were selected for functional screening (Table 1). The pattern of LD was very similar in the 379 samples of European ancestry in phase 1 of the 1000 Genomes project (Additional file 3: Figure S2) [23].

Selection of polymorphisms on the basis of their association with IRF5expression

A recent analysis by our group, in which we combined data from multiple microarray experiments in lymphoblastoid cell lines, highlighted four IRF5 polymorphisms in the best models accounting for IRF5 expression [12]. We selected three of them for the current study. One of them is the CGGGG indel, which is already known to influence IRF5 expression [10, 11]. It was selected as a sort of positive control (Table 1). The second is rs3807306, which had already been selected on the basis of our haplotype analysis. The third is rs17424179, which was not included in our LD analysis because it is 68 kb 3' to IRF5 and we selected it at this stage. The fourth SNP in the microarray-derived model, rs10954213, was not selected because it is in the IRF5 coding sequence and affects mRNA stability.

We selected an additional polymorphism on the basis of IRF5 expression data obtained using the RNA-seq approach (E Alonso-Perez, unpublished data). More details are available in Additional file 1: Note S2. Briefly, the transcriptome of the 60 cell lines of the HapMap CEPH (Utah residents with ancestry from northern and western Europe: CEU) collection was obtained from a published study [24]. After alignment of the RNA-derived sequences, reads corresponding to exons were counted together with the junction reads, and the total was normalized and used as a measure of the expression level. This measure for each cell line was used as the variable to explain the genotypes of 264 IRF5 polymorphisms from the 1000 Genomes project. The most significantly associated polymorphism in this analysis was rs11269962 (P = 0.002), which is a 14-bp indel placed 2.2 kb upstream of the IRF5 gene (Table 1).

Screening of functional polymorphisms with electrophoretic mobility shift assay

We considered it likely that the SNPs of interest could have a functional effect dependent on the differential allele binding of nuclear proteins. EMSA was used as a screening tool for these differences. We used nuclear extracts from a lymphoblastoid cell line in basal conditions. A total of 26 polymorphisms were subjected to this screening (Figure 1): 9 in the rs729302 high-correlation group, 15 from the haplotype selection and 2 on the basis of their association with IRF5 expression (Table 1). In addition, we analyzed the CGGGG indel as a positive control (Additional file 3: Figure S3A).

Only 7 of the 26 polymorphism showed differential allele patterns in EMSA. The first was rs729302 itself, which showed additional binding with the C allele relative to the A allele (Figure 3A). Two others were from the group in tight LD with rs729302. The minor G allele of rs12706860 showed a specific pattern of additional bands with delayed migration when it was compared with the major C allele (Figure 3B). The other SNP of this group, rs13245639, showed specific band slowdown with its minor T allele in comparison with the major C allele (Figure 3C). Also, two of the SNPs selected on the basis of IRF5 protective haplotypes, rs3778754 and rs3807307, showed differential EMSA results. The G allele of rs3778754 showed delayed bands that were much more prominent with the C allele than with the G allele (Figure 3D). Interpretation of this difference was complicated by the lack of specificity of the G allele band, which was attenuated by competition with the unlabelled oligonucleotide but not completely suppressed. In spite of this lack of specificity, other differences between EMSAs with the two alleles led us to select this SNP for further study. The next SNP, rs3807307, showed an additional and specific band shift with the major C allele that migrated faster than the band that was common to the two alleles (Figure 3E). The two polymorphisms selected on the basis of their prominent association with IRF5 expression also showed differential EMSA results. For rs17424179, there were two shifted bands with the two alleles; however, the slowest band was much more intense than the fastest one with the G allele, and this difference was not as marked for the A allele (Figure 3F). For the major allele of rs11269962, the 14-bp insertion, a specific delayed band was observed that was absent with its minor allele (Figure 3G).

These differential patterns of binding with nuclear extracts of the B cell line WIL2 NS were obtained in three or more concordant experiments. None of the other 19 analyzed polymorphisms showed differential results in EMSA. Among them was rs4728142, which has been reported to show differential allele intensity of the delayed band [25]. We analyzed this SNP with particular attention, but we did not detect any difference between the two alleles (Additional file 3: Figure S4).

Luciferase gene reporter assays

To further assess the functionality of the seven polymorphisms that had shown differential EMSA results (Table 1), we performed luciferase gene reporter assays. Each polymorphism was tested with a pair of vectors: one for each allele (148 to 245 bp in length; see Additional file 2: Table S5). They did not contain any other SNP. Inserts for the seven polymorphisms plus those used as positive controls (CGGGG indel) were introduced upstream of the fos minimal promoter in fos-pGL3 vectors. These constructs were transiently transfected into the B-lymphocyte cell line WIL2 NS. The CGGGG indel showed higher expression with the 4x allele than with the 3x allele (Wilcoxon test; P = 0.043), as described previously [11] (Additional file 3: Figure S3B).

Three of the seven newly studied polymorphisms showed differential luciferase expression: (1) rs729302 itself (Figure 4A), (2) rs13245639 (Figure 4C), both from the group of SNPs in tight LD with rs729302, and (3) the 14-bp indel rs11269962 (Figure 4G). The minor C allele of rs729302 showed increased expression relative to the major A allele, but the difference was only borderline (Figure 4A). Constructs with the major C allele of rs13245639 showed a significant increase in luciferase expression compared with the minor T allele (Figure 4C). The construct with the 14-bp indel of rs11269962 showed a significant decrease in luciferase expression relative to the construct with the deletion allele (Figure 4G). These latter results were consistent with the RNA-seq analysis, which showed higher expression with the deletion allele. The remaining four SNPs with differential EMSA results did not show significant differences in the reporter gene assays (Figures 4B, D, E and 4F).

Haplotype assignments of the new functional polymorphisms