Development of the first marmoset-specific DNA microarray (EUMAMA): a new genetic tool for large-scale expression profiling in a non-human primate

Mapping, annotation and submission of sequences to GenBank

A total of 3441 high quality sequences (Phred score: Q15>600; 87.5% of the sequences also passed the threshold of Q20>600 bp; the overall minimal Q20 score was 429 bp) were selected for annotation [13]. After removal of redundant sequences, 3215 high quality unique 3' EST sequences remained which were trimmed to an average length of 678 bases. Annotation of the EST dataset to the RefSeq dataset of human genes allowed unambiguous assignment of gene names for 2250 sequences. A remaining 292 sequences could also be confidently assigned with gene names by mapping to the human genome or cross-species mapping. Of these 2542 (2250 + 292) sequences assigned with a gene name, 610 (19%) sequences could also be annotated with an open reading frame (ORF) (Figure 1A). These ORF-containing sequences were submitted to GenBank [GenBank: EF214838 – EF215447]. The majority of the remaining 673 sequences without an assigned gene name consisted of a group of 642 (20%) which could not be matched to a known mRNA but were mappable on the genome and 31 (1%) which could not mapped to RefSeq genes nor to a mammalian genome (Figure 1A). The latter 31 genes most likely represent novel genes that may be uniquely expressed in the marmoset hippocampus. All 2605 (81%) EST sequences without an ORF were submitted to dbEST [GenBank: EH380242 – EH382846].

Gene selection and construction of the marmoset microarray

The final design of the marmoset microarray consisted of a total of 1649 probe sets of which 97 control probe sets and 1552 probe sets representing 1541 different marmoset transcripts. A total of 71 probe sets were derived from 68 different sequences already available in GenBank, while the remaining 1481 probe sets were derived from 1473 of our proprietary sequences obtained by sequencing marmoset ESTs. The majority (1445) of these 1473 ESTs were derived from the set of 3215 ESTs that have been submitted to GenBank (Figure 1B). Another 28 ESTs had already been selected for representation on the marmoset array based on our initial annotation, but did not pass the criteria we applied for GenBank submission.

The majority of the ESTs represented on the marmoset array assigned with a gene symbol (in total 1251) could be classified according to gene ontology based on biological function (1034) or molecular process (1106), representing a wide variety of different functional classes (Additional file 1).

Expression profiling of hippocampus and several non-neuronal tissues

Correlation coeffcients were calculated based on hybridisation signals on the marmoset microarray for all possible combinations of tissues (Table 3). The highest correlation coefficient (0.96) was obtained when comparing gene expression in both neuronal tissues, i.e. hippocampus and cortex, indicating a high degree of similarity in expression profiles. In contrast, comparing hippocampal gene expression to all peripheral tissues resulted in a severe drop in correlation coefficient (Figure 2), ranging from 0.69 in muscle to 0.50 in liver.

	Hippocampus	Cortex	Fat	Aorta	Ovary	Liver	Muscle	Kidney	Testis
Hippocampus	1	0.95741	0.60714	0.56866	0.54152	0.50381	0.68915	0.54918	0.53917
Cortex		1	0.62146	0.58074	0.54732	0.51323	0.6862	0.74907	0.56225
Fat			1	0.8834	0.822	0.78857	0.7828	0.82475	0.68859
Aorta				1	0.81828	0.76857	0.78129	0.80429	0.67708
Ovary					1	0.76293	0.74486	0.80561	0.70452
Liver						1	0.7282	0.80027	0.65683
Muscle							1	0.74907	0.66837
Kidney								1	0.69856
Testis									1

Half of the probe sets of (767/1552 = 49.4%) were detectable (defined by a present call) in all tissues we tested so far, albeit in most cases at varying levels. Sorting of the expression profiles based on detection call allowed the identification of transcripts expressed exclusively in 1 of the 9 tissues tested. In liver, for example, we identified 7 probe sets with a present or marginal call that had absent calls in the other 8 tissues tested. Three of these 7 probes sets represented cytochrome p-450, enzymes that mediate the oxidative biotransformation of most drugs and are known to be predominantly expressed in liver [14]. Detection of 5 probe sets was restricted to the kidney, one of which represented (prepro)renin. Renin is part of the renin-angiotensin system known to play a major role in volume regulation and is expressed at high levels in the kidney [15]. Sperm protein 17 was exclusively expressed in testis. An explanation for the fact that we find a number of genes that appear to be exclusively expressed in other tissues than brain, while the origin of the ESTs is brain tissue, may be that these genes are predominantly derived from the list of 68 pre-existing GenBank sequences on the array, which are not necessarily expressed in brain. Indeed, this appears to be the case in many instances. For example, the 3 cytochrome p-450 probe sets exclusively expressed in liver are derived from pre-existing GenBank sequences, as well as 4 of the 5 kidney-specific probe sets (including renin) and 3 of the 4 testis-specific probe sets (including sperm protein 17). An additional explanation is that expression levels in brain may be too low to detect using microarrays, but were detectable in our normalised cDNA library of hippocampus.

Detection of a total of 82 probe sets was restricted to neuronal tissue (hippocampus and cortex) (Additional file 3). Among these 82 probe sets were several ion channels (KCNAB3, KCNF1, KCNH3, SCN1a, CACNB1), neurotransmitter receptors (GABRA5, 5-HTR2C, GRIN1), neural cell adhesion molecules (NCAM, OPCML) and diverse components involved in neurotransmitter secretion (SYN2, SNPH), axonogenesis (SLITRK1) and neurotrophic factor signalling (NTRK2). It must be noted that although the sequences have been annotated, the information is not that detailed that for each gene the exact isoform or splice variant is indicated. For example, a closer look at the data of CACNB1 indicates that only isoform 1 is targeted by our sequence, which indeed is expressed in brain, while isoform 2 is expressed in muscle. The users of the marmoset microarray have access to the exact probe set sequences as well as the EST sequences, allowing unambiguous determination of which part of the transcript is targeted and which splice variants or isoforms are detected.

Ninety probe sets were not detected in any of the tissues tested, which may be an intrinsic property of the selected probes or due to potential sequencing errors resulting in mismatches of the probe sets with the labeled mRNA (Additional file 3).

Generation of marmoset cDNA sequences

Total RNA was isolated from hippocampi of 5 adult marmosets, both males and females ranging from 3 to 10 years in age, After checking the quality of the isolated RNA it was pooled and used to construct a custom normalised cDNA library in pCMVSport6.1 (Invitrogen, Carlsbad, USA). For the development of the marmoset microarray a total of 23,300 clones were gridded in 61 384-well plates. In total approximately 5700 sequencing reactions were performed at Baseclear (Leiden, The Netherlands) using an anchored oligo-dT primer on 15 randomly selected 384-well plates to generate partial 3' EST sequences.

RNA preparation and labelling

Total RNA was isolated using TRIzol^® reagent (Invitrogen Life Technologies, Carlsbad, CA, USA) according to the manufacturer's instructions. After isolation, total RNA was purified using the QIAGEN RNEasy^® Mini Kit RNA Cleanup procedure (QIAGEN Inc., Valencia, CA, USA). RNA quality was assessed with the LabChip^® RNA 6000 Nano Assay on the 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA, USA). Per RNA sample, 5 μg was used as input into the MessageAmp™ II-Biotin Enhanced Single Round aRNA Amplification Kit (AM1791, Ambion, Austin, USA) according to the manufacturer's instructions. Briefly, total RNA was converted to double-stranded cDNA followed by in vitro transcription to synthesize biotin-labeled aRNA. After purification the size distribution of the aRNA was checked with the LabChip^® RNA 6000 Nano Assay on the 2100 Bioanalyzer (Agilent Technologies) and the quantity was assessed by small sample spectophotometry on a Nanodrop ND-1000 Spectrophotometer (NanoDrop Technologies, Wilmington, Delaware USA).

Hybridisations

The biotinylated aRNA was fragmented and 15 μg was hybridised to marmoset arrays at the Leiden Genome Technology Center (LGTC^®), Leiden University, The Netherlands, according to the manufacturer's recommendations (Affymetrix, Santa Clara, USA).

Data analysis

Affymetrix GeneChip^® Operating Software (GCOS) was used to estimate signal intensities plus signal reliabilities. Data was median normalized in excel and subsequently scatterplots were constructed and correlation coefficients calculated using R version 2.2.0 [33]. Probe sets assigned with a gene symbol were functionally classified according to biological process and molecular function using Netaffx Analysis Center [34].

Mapping and annotation of EST sequences

In order to obtain high-quality data the EST sequences were trimmed before further analysis. At the 3' end, at least 15% of each sequence was removed so that final DNA fragment lengths were 700 bases. In addition, 20 bases from the 5' end of the EST sequences were trimmed. For the annotation of the EST dataset, the sequences were reverse complemented and compared to the RefSeq dataset of human genes (NCBI Refseq build 36.1) using a locally installed version of the Blat software [35]. Remaining sequences that could not be unambiguously assigned to a gene name were mapped to the human genome sequence (NCBI Refseq assembly version 36.1) according to the standard parameters of the software. ESTs with a hit on the genome were annotated with the gene name of the gene which is located at the respective chromosome region. Ambiguities, e.g. if an EST sequence covered more than one gene, were resolved manually. Finally, a Blat search with increasing sensitivity parameters was performed to map ESTs across species.

The annotation procedure described above allowed a reliable annotation of EST data even if sequencing errors would cause frame-shifts in coding regions. In order to annotate open reading frames (ORFs) associated with the ESTs, ORFs with a length of at least 100 bases starting from the 5' end of the reverse complementary sequences were extracted. The translated amino acid sequences were mapped against the human protein sequences from UniProt (Release 50 of May 2006, [36]) using Blat [35]. All sequences without any hit on a protein sequence starting at the first position ofthe 5' end where rejected. The remaining candidate ORFs were first validated by a comparison of the gene names from the UniProt hits to the gene names of the above determined RefSeq hits. To be included as bona fide ORFs, the translated ORF sequences had to fulfil an additional criterion. The length of the peptide sequences derived from the ESTs had to cover at least 80% of the peptide sequence length which can be expected from the corresponding part of the protein sequence from UniProt.