Development of a chicken 5 K microarray targeted towards immune function

Construction of the array

Genes on the array

Analysis of the immune clones

All the sequences of the clones on the array were subject to Blast homology searches against the SwissProt and TREMBL databases using a cut-off value of 1e^-10. Using this means of detection, many known immune-related molecules were identified, including cytokines, interferons, interleukins, transcription factors, receptors, cell differentiation antigens, MHC molecules and genes for proteins belonging to the TOLL receptor pathway. Proteins homologous to hypothetical human proteins and mouse cDNAs were also identified.

Sequences, which gave no Blast homology to anything in the nucleotide or protein databases, accounted for about 38% of the clones. Either the search parameters were too stringent to identify these genes or the chicken sequence was sufficiently divergent to be undetectable in a standard Blast search. This is a common feature of immune-related genes, and it is often very difficult to identify such genes by sequence homology to mammalian homologues. Some of these sequences may also represent non-conserved 3' UTR regions of genes. This set of clones may also include genes that have never been identified before and are not represented in the sequence databases. Further, more detailed analysis of these sequences can sometimes help elucidate the nature of the gene in question. Protein sequences can be predicted from the EST nucleotide sequence using programs such as ESTscan [[16] and [17]], which takes in to account sequencing errors and thus potential frame-shift mutations which are often present when there is only one EST sequence available for study. Conserved motifs and domains can then sometimes be identified for example, using the Pfam database [18], which is a large collection of multiple sequence alignments and hidden Markov models covering many common protein domains and families. PSI-Blast searches can also help identify to which type of family a gene will belong.

Gene ontology (GO) annotations

In order to try and elucidate the function of the genes on the array further, we tried to assign as much annotation to the sequences as possible. GO annotations were assigned to some sequences after searching the GGI and UMIST databases [19], while other annotation was derived from hits to orthologous human sequences from the ENSEMBL [20] and GENSCAN [21] databases, as described in the 'methods' section. Having annotation derived from orthologous human genes means that cross-species comparisons between chicken and human array data may be possible. A search of the ENSEMBL database provided information on 2,292 GO-term associations, the GGGI database 1,542 and GENSCAN 566, while the UMIST full-length cDNA database provided a further 365 annotations. The sequences on the array cover a total of 227 GO terms, with 73% of all the sequences having at least one GO entry assigned to it. The available annotation for the array sequences is broken down as follows: 52% of genes have a 'cellular component' term assigned, 60% have 'molecular function' and 56% of sequences have the 'biological process' described. 83% of all the genes on the array have some kind of gene description and after searching each sequence against the sequences in the Ensembl chicken genome collection (July 2005 genebuild [22]), 78% of sequences were found to have a known chromosomal location. Now that all these sequences have been added to GenBank and thus have an accession number which can be directly linked into the ENSEMBL databases (work currently underway), obtaining comprehensive, up-to-date annotation data will become much easier.

A file showing the complete annotation for all the sequences on the array is available as supplementary material (Additional file 1). However, Additional file 2 provides an overview of the broad functional classes that are represented by the genes on the array. These are based on more general GO annotations derived from the GO-slims database at EBI, and allow us an insight into the different classes of genes present on the array without having to look at detailed functional annotation for each individual gene.

Annotation is also available for some (9,137) of the ESTs in the UMIST collection. By comparing the relevant GO slims [23] terms for the sequences in this collection with those present on our array, we are able to see which types of genes appear to be enriched in our set, compared with a larger, more general collection of EST sequences. As can be seen (shown in bold) in Additional file 2, certain classes of gene appear to be more highly represented. For instance, genes involved in protein transport are more abundant in our set of clones, as are those involved in the response to stimulus. This is consistent with our attempts to pre-select for higher numbers of genes involved in the immune system.

Quality of the array

To assess the quality of the array, various hybridization comparisons were undertaken. Three different conditions were addressed: 1). self v self 2). biological replicate A v biological replicate B and 3). Control sample v activated sample. Dye swap experiments were also carried out for conditions 2 and 3. The 'self' sample was a reference RNA consisting of a pool of various chicken lung samples. The biological replicates were lung samples from two 6-week-old chickens that had not been treated or challenged in any way. In the third group of hybridisations, the 'control' sample was from a similarly, untreated bird and the 'activated' sample was obtained from the lungs of a bird that had been challenged with the avian influenza strain H9N2 five days previously. The graphs in Fig 1 show the tight correlation between self/self (R² = 0.9273) and between replicates (R² = 0.8766), whereas a much higher level of variance is seen when an activated sample is compared against a control (R² = 0.7601).

The boxplots in Fig 2 also demonstrate the differing variances between the comparisons. The greatest variance is shown for the activated animals compared with the controls as would be expected. Regression analysis for each of the data sets confirm the increased variance with correlation coefficients of r = 0.872 for activated samples, r = 0.936 for replicate samples and r = 0.963 for self/self sample data sets.

Using the array

This array is available from the Ark-Genomics resource facility at Roslin Institute, providing an immune-focused array which, for anyone interested in immune-research, offers a much more cost-effective and time-saving platform for gene expression experiments, instead of using the large oligo arrays which have thousands more genes, many of which will be of no interest. Analysis of data is also thus much easier and far less time-consuming. Information on the array has been deposited in ArrayExpress (Accession: A-MEXP-307) [[24] and [25]] (Additional file 1) and very soon all the sequences will be submitted to the Ensembl database with links to all the GO annotation information in the GOA database [26].

Sample collection

Eight groups of 38 chickens (3-week-old) were vaccinated with two different vaccine regimes. The eight groups were males and females of a commercial line of hybrid broiler (Ross 306, Aviagen, Newbridge, Midlothian, UK) and layer (Lohman Brown, Lohmann Tierzucht, Cuxhaven Germany) chicks given one of the two vaccination schemes. Group 1 were given vaccines for E. coli (0.5 ml in left breast muscle), ND and IBDV (0.5 ml in right breast muscle) formulated in alum-gel and oil-based immuno-potentiators. Intramuscular injections were given to ensure that all the birds were given an equal dose. Group 2 vaccines consisted of Paracox 8 [Eimeria sp.] (0.1 ml in drinking water), Nobilis Rismavac-CA126 [MD] (0.2 ml intramuscularly in leg) and Salenovac [S. enteritidis] (0.5 ml intramuscularly in leg). Tissue samples were obtained (unvaccinated); 5 hr, 24 hr, 72 hr and 7 days post vaccination. Samples from groups of 5 birds were pooled. Tissues collected were Bursa, spleen, Peyers patch and thymus. Tissue from Bursa, spleen and Peyers patch were pooled to make the 'B-cell' libraries and the thymus tissue was used to construct the 'T-cell' libraries. The tissues and time points chosen were in order to try and maximise the number of immune-related transcripts, including those which may only be expressed transiently. All experimental protocols were authorized under the UK Animals (Scientific Procedures) Act, 1986.

Library construction

Six libraries were constructed at Incyte Genomics (Palo Alto, CA): a standard and 2 normalized Bursa/spleen/Peyers patch libraries and a standard and 2 normalized thymus libraries. cDNA synthesis was initiated using an oligo (dT) primer, using methylated C in the first strand synthesis reaction. Following this first strand reaction, double-stranded cDNA was blunted, ligated to NotI adapters, digested with EcoRI, size-selected, and cloned into the NotI and EcoRI compatible sites of a custom modified MCS of the pBluescript (KS+) vector. Normalization was done in two rounds using conditions adapted from [[27] and [28]] except that a significantly longer re-annealing hybridization was used. Around 10,000 clones were then sequenced at the Sanger Institute according to their protocols. Using the T7 primer, sequence was generated from the 5' end of each clone by the dideoxy chain termination method using an ABI 3700 sequence analyser (Applied Biosystems, Foster City, CA).

EST sequence analysis

Bioinformatic analysis commenced with 10,173 sequences. After eliminating poor quality sequence and repeats, 9,434 of these sequences remained after screening with phred [29], RepeatMasker [30], Crossmatch [31] and XNUN [32]. Certain unwanted sequences were then identified after using the Blast algorithm [[33] and [34]] and screening the results for specific keywords. These included 'ribosomal', 'mitochondrial', 'Newcastle', 'Mareks', 'Eimeria', 'Salmonella' and 'E. coli'. 8,154 sequences passed these criteria. These sequences were then clustered against the existing UMIST and EMBL chicken EST sequences using TIGR's clustering tool, tgicl [35]. This resulted in 3,845 clusters which contained one or more sequence from our libraries and 1,959 singletons. The following clones were chosen for inclusion on the array: 3,770 cluster representatives, 1,067 singletons and 157 reference immune genes: 93 clones from the UMIST collection, 41 from our immune libraries, 21 clones from the Delaware set [36] and 2 clones courtesy of R. Zoorob (CNRS, France) (Table 2).

Construction of the array

The immune array was constructed from 4994 chicken EST clones plus 196 control elements (landing lights (positional controls), GAPDH, gamma actin, salmon sperm DNA, calf thymus DNA, chicken and bovine genomic DNA and a variety of spotting buffers). Plasmid DNA was prepared using MagAttract 96 Miniprep chemistry on a Biorobot 8000 platform (Qiagen Ltd., Crawley, UK), and the cDNA inserts were amplified using CGATTAAGTTGGGTAACGC (fwd) and CAATTTCACACAGGAAACAG (rev) in 50 ul reactions using 1 ul of DNA as a template. Amplified DNA was purified by Multiscreen 384 well PCR purification plates (Millipore, Watford, UK) on a Multiprobe II liquid handling platform (Perkin Elmer, Beaconsfield, UK) and the reactions confirmed by agarose gel electrophoresis and quantified by Picogreen assay (Molecular Probes, Invitrogen, Paisley, UK) on a Flouroskan Ascent flourescent plate reader (Thermo Life Science, Basingstoke, UK). DNA was resuspended to 150 ng/ul in spot buffer (150 mM Sodium phosphate, 0.01% SDS) before being spotted in duplicate on to amino-silane coated slides (CMT-GAPSII, Corning, Schiphol-Rijk, The Netherlands) using a Biorobotics MicroGrid II spotter (Genomic Solutions, Huntingdon, UK). Slides were then treated using succinic anhydride and 1-methyl-2-pyrrolidinone (Sigma, Poole, UK) to block unbound amino groups, followed by a wash in 95°C MilliQ water before hybridisation.

RNA preparation and labelling

Total RNA was isolated from lung tissue using a Trizol extraction according to the manufacturer's protocol (Invitrogen, Paisley, UK) and subsequently purified using the RNeasy Midi RNA Purification kit (Qiagen Ltd., Crawley, UK). RNA concentration was determined spectrophotometrically and RNA quality was determined using an Agilent 2100 Bioanalyser (Agilent Technologies, Waldbronn, Germany). Cy3 or Cy5 was incorporated into each sample using the Fairplay labelling kit (Stratagene, La Jolla, CA) and the labelled cDNA cleaned-up after passage through DyeEx columns (Qiagen Ltd., Crawley, UK). Labelling efficiency was determined by running 0.5 μl of each sample on a 1% agarose gel and measuring the intensity of fluorescence on a GeneTac LS IV scanner (Genomic Solutions, Huntingdon, UK).

Hybridizations

Microarray hybridizations were carried out overnight using a GeneTAC automated hybridization system [37] (Genomic Solutions, Huntingdon, UK). Hybridizations (125 μl) were carried out in Genomic Solutions hybridization solution (Cat. no. RP#0025) in a stepped hybridization: 55°C for 3 hr, 50°C for 3 hr and then 45°C for 12 hr. Slides were then washed in Genomic Solutions wash buffers (Cat. nos. CS#0038, CS#0039 and CS#0040). Upon removal from the hybridization stations, slides were washed for 1 min in Post-Wash buffer (CS#0040) and a further minute in isopropanol, followed by centrifugation at 1000 rpm for 6 min. Dried slides were scanned in a Scanarray 5000 scanner (GSI Lumonics, Rugby, UK) fitted with Cy3 and Cy5 filters.

Data analysis

To indicate the suitability of the new array to discriminate the differences in the experimental treatments, hybridizations comparing samples with controls and controls with controls were performed. Control (vehicle treated) animals were compared with immunologically challenged animals (activated slides) and control animals were also compared with other control individuals (replicate slides). The same animal was also compared with itself (self/self). Each comparison was completed in duplicate and with a dye flip. Dye-swaps are carried out in order to deal with any residual dye-bias remaining after labelling. However, this is generally not a problem, due to the indirect labelling method employed. Data was extracted from the slide using Bluefuse software (BlueGnome, Cambridge, UK). Features with poor confidence information (confidence <0.30, flagged D and E) were eliminated from the analysis. M v A plots [where M = log₂ (Cy5/Cy3) and A = 1/2*(log₂(Cy5) + log₂(Cy3)] of the data for each slide (data not shown) were suitably linear to require only a simple global normalisation of the data. Data from slides of similar treatments was pooled and a boxplot produced for each comparison (Genstat v8.1, VSN International Ltd., Hemel Hempstead, Herts, UK).

Databases and sequence sources

Ensembl and Genscan predicted genes/peptide sequences for the chicken genome assembly (March 2004) were downloaded from the Ensembl database using Ensmart or the UCSC table browser [38]. Chicken EST sequences were downloaded from the TIGR Gallus gallus gene index (GGGI) [release 10.0] [[39] and [40]]. Chicken full-length cDNA sequences were downloaded from the UMIST www site (Sept 2004). Ensembl predicted peptide sequences for the human genome assembly (May 2004) were downloaded from the Ensembl database using Ensmart or the UCSC table browser.

Mapping array probes to chicken ESTs, cDNAs, genes and genome

Unique ESTs used to create the immune array were mapped to chicken cDNAs, ESTs, genes or the chicken genome assembly using NCBI Blastn (version 2.2.11). Identity was defined with > 95% sequence identity over 100-bp and then taking the top-scoring match to each EST to provide a unique sequence assignment. All repeats and low-complexity sequences were masked using RepeatMasker (version 3.1.0).

Definition of Gene Ontology terms and Gene Descriptions for array probes

Gene Ontology (GO) annotations [41] were all based on database hits in sequence similarity searches using Blastn. GO annotations were automatically transferred from these database records to the array probe entries. GO annotations were available for GGGI and UMIST EST/cDNA sequences. For chicken Ensembl or Genscan gene predictions, GO annotations were based on orthologous human peptide sequences. Orthologues were defined based on two cycles of Blastp between human and chicken proteins. An E_value cut off of less than 10^-4, with the subject and query databases swapped between runs. By comparing E_values mutually best proteins pairs were selected as orthologues. When E_values were equal, bits score and sequence coverage were used as tiebreakers to select the top hit. For each array probe associated GO terms and a unique gene description was transferred from the orthologous database record. Finally a Perl script was used to create a non-redundant set of probe to GO records.

Frequency of GO and GO-Slim terms

GO terms (version 3.2.16) were downloaded from the Gene Ontology www site. More general GO terms were assigned using GoaSlim_map (June 2005) available from the GOA www site at EBI. The GO-Slim terms allowed us to estimate e.g. the frequency of array probes associated with the biological process Metabolism (GO:0008152).

Data processing

Perl scripts (version 5.8.5) and SQL were used throughout to manipulate and filter data sets.