text:
miRNAminer: A tool for homologous microRNA gene search

Background MicroRNAs (miRNAs), present in most metazoans, are small non-coding RNAs that control gene expression by negatively regulating translation through binding to the 3'UTR of mRNA transcripts. Previously, experimental and computational methods were used to construct miRNA gene repositories agreeing with careful submission guidelines. Results An algorithm we developed – miRNAminer – is used for homologous conserved miRNA gene search in several animal species. Given a search query, candidate homologs from different species are tested for their known miRNA properties, such as secondary structure, energy and alignment and conservation, in order to asses their fidelity. When applying miRNAminer on seven mammalian species we identified several hundreds of high-confidence homologous miRNAs increasing the total collection of (miRbase) miRNAs, in these species, by more than 50%. miRNAminer uses stringent criteria and exhibits high sensitivity and specificity. Conclusion We present – miRNAminer – the first web-server for homologous miRNA gene search in animals. miRNAminer can be used to identify conserved homolog miRNA genes and can also be used prior to depositing miRNAs in public databases. miRNAminer is available at .

Background
MicroRNAs (miRNAs) are short, ~22 nt non-coding RNAs that control gene expression. miRNAs bind to the 3'UTRs of their regulated mRNA transcripts to facilitate mRNA degradation or translation inhibition [1,2]. miRNAs are present in most metazoans and are thought to regulate a diverse range of biological processes [3,4]. miRNA genes' evolution is unique since they mostly emerge through duplication events [5]; exhibit most often unidirectional evolution [6]; are generally gained and not lost [7]; show several cases of rapid evolution in primates [8,9]; are rarely changed due to functional constraints [10]; and, show relatively rare evolutionary acquisition events (accounted by their relatively small numbers). miRNA predictions based on computational methods, which followed initial extensive cloning efforts, are based on the secondary structure of the miRNA, its phylogenetic conservation and thermodynamic stability [11,12]. miRNA gene repositories are constantly expanding giving rise to more than 3500 reported miRNAs in more than 30 animal species (Sanger miRbase database, Version 10.0 [13,14]). However, it is apparent that even this comprehensive repository is far from completion, accounted by the very few miRNAs listed for dog (6) and chimpanzee (83), compared to human (533), to name two examples. Since these differences cannot be accounted merely by species-specific miRNAs, we saw the need for a computational tool for miRNA homologous searches.

Implementation
We present miRNAminer, a tool for automatic identification of homolog miRNAs based on a given user defined query miRNA. The tool exploits numerous characteristics of miRNAs: high conservation of precursor sequences, very high conservation of mature sequences (particularly in the seed region, nt 2-8 [15]), and hairpin secondary structure with high folding energy and base pairing. miR-NAminer first uses BLAST [16] to select candidate matches and ranks them according to their e-values. Then it employs a series of rigorous filters to improve specificity.
An input query consists of a precursor miRNA, mature miRNA, a set of filter threshold values and the number of best-fitted results requested in the output. We designed miRNAminer's algorithm to maximize specificity of matches. This is because the designed application of miR-NAminer is to identify homolog matches after a miRNA has been experimentally confirmed. We estimated the default values presented below so that each filter by itself selects 95% of known miRNAs in training genomes (criteria was also based on [17]). miRNAminer's algorithm follows these steps: (i) Use BLAST [16] to find matches in target genomes (the whole precursor miRNA from the query is used); (ii) Filter with e-value threshold (default 0.05 per chromosome); (iii) Extend the match by adding flanking nucleotides (default 50) up-and down-stream from the match (Ensembl genome database; [18]). Examine all possible extensions of the match within threshold length (default min 70 nt, max 180 nt); (iv) Filter with RNA secondary folding energy threshold (default -25 kcal/mole; RNAfold with options "-p -d2 -noLP" [19]); (v) Filter with minimal base-pairing threshold (default 55% pairing; with 20 gap penalty and 0.5 extension penalty); (vi) Filter with requirement for hairpin-shape secondary structure; (vii) Filter with alignment of precursor sequences (default 56% identity); (viii) Filter with alignment of mature miRNA sequences (default 80% identity); (ix) Filter with maximum number of mismatches in mature miRNA sequences (default 3 nt); (x) Filter with conservation of seed (2-8 nt, required 100% conservation [15]); (xi) Filter with position of mature miRNA on the hairpin (max 4 nt overlap of mature sequence and hairpin loop). miRNAminer's output includes detailed analysis of the identified genomic region(s) that passed the selected threshold criteria ( Figure 1). Currently, miRNAminer supports searches in 10 metazoan genomes. We will regularly add additional genomes upon their release. After the query is issued, results are usually available within a minute (though this depends on the number of results requested) and can either be viewed on the screen or requested to be sent by email.

Results
We used miRNAminer to perform a comprehensive homology search for miRNA precursors in seven species (human, chimpanzee, mouse, rat, dog, cow and opossum). For the search, we used all 2925 vertebrate miRNAs listed in the Sanger miRNA registry (release 9.0 of October 2006). Table 1 shows the summary information of miR-NAs listed in the Sanger registry and of new, or non-registered, miRNAs identified by our method. We identified 790 non-miRbase registered miRNAs with major contributions to chimp (P. troglodytes), dog (C. familiaris) and cow (B. taurus), vastly increasing their known miRNA repertoire and possibly opening new research facets in these species (see Additional file 1[20]). Table 2 presents the miRNA candidates that our method identifies in human (H. sapiens). It is of interest that 22 new candidate miRNAs in human were identified, despite many previous exhaustive human miRNA identification studies [8,21], possibly due to recent non-human/primate miRNAs identification [22], updated assembly of the human genome (Ensembl H. sapiens genome release 42 used here was updated in October 2006), and the modification of search parameters as implemented in miRNAminer.
We compared our predicted miRNAs (in the human species; Table 2) to other prediction methods. We found that 18 and 36% of our miRNAs are contained within RNAmicro [23,5] and Berezikov 2005 [24] databases, respectively, out of almost 3000 and 1000 miRNAs in each set, respectively. The overlap is not extensive however even when miRNAs derived from algorithms using very similar search parameters are compared only about 50% overlap is seen [25]. One of our identified miRNAs, which is not reported by any other study, was recently identified experimentally ( [26] present in miRbase version 10.0 and not 9.0) increasing confidence in our unique miRNAs. Notably, even though miRNAs from all species were used to search for human homologs, the candidate miRNAs discovered are homologs to genes in two species only, M. musculus and R. norvegicus, indicating better miRbase coverage for mouse and rat than for other species. Two examples of the non-registered human homolog miRNA genes are presented in Figure 2. miRbase mouse miR-764 sequence, which has no known registered homologs was used as input for miRNAminer search (with default parameters). The output reported a homolog (presumably hsa-miR-764; see Table 2), which is located in the second intron of Human serotonin receptor 2C (HTR2C; NM 000868; Figure 2A). The mouse miRNA homolog is located in an intron of the same gene (HTR2C; NM 008312) suggesting an evolutionary conserved co-expression of miRNA and its host gene [27][28][29]. Second in the list of non-registered human homolog miRNA genes (sorted according to BLAST e-values) is miR-763. This miRNA spans the longest complementary sequence out of the list, has the lowest RNA folding energy (Table 2) and shows high conservation between many species ( Figure  2C). Interestingly, human miR-763 is harboured in an intron of the high mobility group AT-hook 2 oncogene (HMGA2; for review see [30,31]). Recent disrupted interplay between miRNAs and HMGA2 showed an increase in oncogenesis [32][33][34]. To regulate their targets, miRNAs bind to 'seed' regions in the 3'UTR, typically 6-7 nt long (nt 2 to nt 7 or 8 of the miRNA [15,35,36]; also see [37]). miR-763, possibly also co-expressed with its host gene [27][28][29], has a conserved binding site for its own harboured miRNA (nt 2 to 8 of the miRNA binds position 2192 in the 3'UTR which is conserved in human/mouse/ rat). It is tempting to speculate a negative feedback regulatory role of newly identified human miR-763 and its oncogenic host when co-expressed in the same spatiotemporal context. To this end data from Expressed A sample output of miRNAminer Figure 1 A sample output of miRNAminer. Indicated are: search start and end times (rows 1 and 24, respectively) and the species and assembly searched (rows 2-3); whether a match (a miRNA homolog), or matches, that passed the input criteria were found (rows 4-5); information about the quality of the homolog miRNA match such as BLAST e-value, genomic location, sequence, RNA fold and energy, pairing, length and alignment with input sequence (rows [6][7][8][9][10][11][12][13][14][15][16][17][18][19]. A hyperlink to the genomic locus of the miRNA homolog is also provided through Ensembl ContigView [18] or UCSC Genome Browser [43] (rows 20-23) and a copy of the users' input data (row 25). The column 'miRNA registry 9.0' shows the number of miRNAs listed for the given species in the Sanger miRNA registry release 9.0. The column 'Newly identified' shows the number of miRNAs candidates identified by miRNAminer.
Sequence Tags (ESTs) supports this possibility ( Figure  2C). Other identified miRNAs presented in Table 2 show high species conservation (for example, miR-670) or are located in exons (for example, miR-711) or exon-intron junctions (for example, miR-762). Interestingly, in a recent study involving deep sequencing, four of our human predicted miRNAs were confirmed (miR-760, 708, 543, and 665 [38], available in miRbase version 10.0). To conclusively confirm the presence of the identified candidates in the studied species, an experimental verification is required. However, the candidates identified by our method are close homologs to known miRNAs and as such are not required to meet as stringent criteria to be annotated as novel miRNAs [17]. In this study we looked at homolog genes which are genes related to each other by descent from a common ancestral DNA sequence. We do not segregate between orthologs, genes in different species that evolved from a common ancestral gene by speciation, and paralogs, genes separated by the event of genetic duplication. We cannot also rule out that similar miRNAs in different species have developed independently [39]. Our tool, which is based on evolutionary conservation, can only detect evolutionarily conserved miRNA genes. We are currently improving our algorithm to include multiple alignments of vertebrate miRNA sequences in order to better refine the boundaries of the miRNA precursor sequence.  The results are sorted according to e-values (based on BLAST search; [16]). Lengths refer to lengths of miRNA precursors. For each result, the table shows position, RNA fold quality and the identity ratio of the discovered miRNA candidates to the homolog in the original species. The existence of two occurrences for miR-322 and miR-543 (italicized) indicates that two origins (mouse and rat) of the same miRNA were found in Homo sapiens.
additional miRNA in human and none in mouse. This, however, might change when run in combination with other modified parameters. Altogether each of the modified parameters listed above will result, independently, in an average miRNA increase of 11% when tested on seven mammalian species.
We estimated miRNAminer's sensitivity (Table 3) and specificity. The sensitivity, on seven mammalian species, is 0.88. Sensitivity for a species is the portion of the species' miRNAs with known homologs that are detected by miRNAminer using miRNAs from all other species. We used only miRNAs which miRbase lists for more than one species. Sensitivity measures are higher in chimp (0.94), mouse (0.88) and rat (0.91) than in human (0.85). To Two examples of a non-miRbase registered miRNA identified using our miRNAminer web-server Figure 2 Two examples of a non-miRbase registered miRNA identified using our miRNAminer web-server. (A) Human miR-764 was identified using miRbase mouse miR-764 sequence as input (and default parameters) for miRNAminer search. The output reported a homolog (presumably hsa-miR-764), which is located in the second intron of human serotonin receptor 2C (HTR2C; NM 000868). The mouse miRNA homolog is located in an intron of the same gene (HTR2C; NM 008312) suggesting an evolutionary conserved co-expression of miRNA and its host gene [27][28][29]. High conservation is seen in this region (mountain-like graph derived from UCSC Genome Browser 17 species multiZ alignment; [43]). Black rectangles represent exons (shorter rectangles in C are UTRs), lines are introns and dark-grey rectangles are miRNA genes. (B) RNA secondary structure of both the identified human (top) and mouse (bottom) miR-764 exhibit similar thermodynamic stability (41.8/49.9 kcal/mol, respectively) and structures (mature miRNA region is underlined). Human miR-764 homolog was also identified by Berezikov [21]. (C) Non-registered (miRbase) human miR-763 is highly conserved among vertebrate species and can potentially bind its own host gene. On top; a schematic non-scaled representation of the HMGA2 transcript (NM 003483; human miR-763 is in dark-grey; conservation plot as shown in A). Expressed Sequence Tags (ESTs; light-grey bars) are evidence for the expression of this particular genomic region. ESTs from top to bottom: BM715067 (isolated from eye-related tissue); BJ997562 (isolated from wilms tumor tissue); BU39975 (isolated from eye-related tissue); AI935081 (tissue source unknown). On the right; the potential binding site of miR-763 in HMGA2 3'UTR (nt 2-8 of the miRNA; positions 2192-2198) is conserved to human, mouse and rat.
estimate specificity, we used miRNAminer to search for miRNA homologs in C. elegans, which has a large evolutionary distance from the studied mammals. We treated as false positives all hits reported by miRNAminer that were not identified as homologs by previous studies. This conservative treatment may over-approximate the number of false positives. Using 1375 miRNAs from the seven studied mammalian species, miRNAminer detected, in C. elegans, two known homologs (let-7 and mir-124) and reported only five false positives.

Conclusion
Several approaches to identify miRNA homologs have been previously described, both in plants [40], and in animals [5,41,42]. However, the only tool that is available as a web service, microHARVESTER [40], is targeted for plants. miRNAminer is the first available miRNA gene homolog search tool for animal genomes.

Additional file 1
Non-miRbase registered miRNAs. A list of 790 miRNAs that were identified using miRNAminer. These miRNAs add more than 50% to the total count of miRNAs in the seven mammalian species tested: human, chimpanzee, mouse, rat, dog, cow and opossum, and are available at: http:// web.mit.edu/nshomron/www/miRNAminer_SM1.zip Click here for file [http://www.biomedcentral.com/content/supplementary/1471-2105-9-39-S1.zip] The second (third) column presents the number of miRNAs with homologs that miRNAminer detects (does not detect) in each species.


extfieldsofstudy: Biology Medicine

s2fieldsofstudy: Biology Medicine

text:
GPRC6A Null Mice Exhibit Osteopenia, Feminization and Metabolic Syndrome

Background GPRC6A is a widely expressed orphan G-protein coupled receptor that senses extracellular amino acids, osteocalcin and divalent cations in vitro. The physiological functions of GPRC6A are unknown. Methods/Principal Findings In this study, we created and characterized the phenotype of GPRC6A −/− mice. We observed complex metabolic abnormalities in GPRC6A −/− mice involving multiple organ systems that express GPRC6A, including bone, kidney, testes, and liver. GPRC6A −/− mice exhibited hepatic steatosis, hyperglycemia, glucose intolerance, and insulin resistance. In addition, we observed high expression of GPRC6A in Leydig cells in the testis. Ablation of GPRC6A resulted in feminization of male GPRC6A −/− mice in association with decreased lean body mass, increased fat mass, increased circulating levels of estradiol, and reduced levels of testosterone. GPRC6A was also highly expressed in kidney proximal and distal tubules, and GPRC6A−/− mice exhibited increments in urine Ca/Cr and PO4/Cr ratios as well as low molecular weight proteinuria. Finally, GPRC6A −/− mice exhibited a decrease in bone mineral density (BMD) in association with impaired mineralization of bone. Conclusions/Significance GPRC6A−/− mice have a metabolic syndrome characterized by defective osteoblast-mediated bone mineralization, abnormal renal handling of calcium and phosphorus, fatty liver, glucose intolerance and disordered steroidogenesis. These findings suggest the overall function of GPRC6A may be to coordinate the anabolic responses of multiple tissues through the sensing of extracellular amino acids, osteocalcin and divalent cations.

Introduction
GPRC6A is a recently identified member of family C of G protein-coupled receptors (GPCRs) that is most closely related to the calcium-sensing receptor CASR [1,2,3]. Structural homologies and conservation of specific domains in members of this family of receptors suggest an evolutionary link between extracellular calcium and amino acid-sensing [4,5,6]. Indeed, GPRC6A has recently been shown to sense extracellular cations and amino acids and to require both extracellular cations and amino acids for optimal stimulation in vitro [3]. This dual sensitivity of GPRC6A to both divalent cations and amino acids is analogous to the related receptor CASR [7]. Compared to CASR, much higher extracellular calcium concentrations are needed to activate GPRC6A [3,8], and some studies suggest that cations may only be allosteric modulators of GPRC6A [9], whereas others show cationdependent activation of GPRC6A [3]. The calcimimetic NPS-R578, an allosteric modulator of CASR [10,11,12], and osteocalcin, a bone derived calcium binding protein, both enhance the functional responses of GPRC6A to extracellular calcium in vitro [3]. The physiologically relevant ligands for and biological function of GPRC6A remain to be determined [13].
GPRC6A is broadly expressed in many tissues and organs, including lung, liver, spleen, heart, kidney, skeletal muscle, testis, brain and bone [1,2,3,14]. The amino acid, osteocalcin, and divalent calcium ligand specificity of this receptor and its wide tissue distribution implicate GPRC6A multiple processes. For example, GPRC6A may be a candidate for the elusive extracellular calcium-sensing mechanism known to be present in osteoblasts [3,15,16,17,18], which respond to high local Ca 2+ concentrations (in the range of 8 to 40 mM) [19], amino acids and osteocalcin in the bone microenvironment. GPRC6A is also a candidate for the putative osteocalcin receptor regulating energy metabolism [20]. To determine the function of GPRC6A, we created mice lacking this receptor. We found a complex multiorgan, metabolic-like syndrome in GPRC6A 2/2 mice that suggests that GPRC6A is involved in nutritional pathways coordinating the metabolic activity of multiple tissues in response to changes in extracellular amino acids and divalent cations.

Creation and Gross Phenotype of GPRC6A Deficient Mice
We selectively deleted exon 2 of the mouse GPRC6A gene (Fig. 1A). Wild-type GPRC6A +/+ , heterozygous GPRC6A +/2 , and homozygous GPRC6A 2/2 mice were genotyped by PCR (Fig. 1B) and each genotype was found to be born at the expected Mendelian frequencies from heterozygous breeding pairs, however, breeding pairs consisting of GPRC6A 2/2 6GPRC6A 2/2 had litter sizes that were roughly half that of GPRC6A 2/2 6GPRC6A +/+ or GPRC6A +/2 and GPRC6A +/2 breeding pairs, respectively (Data not shown). In addition, full-length GPRC6A transcripts and proteins were documented to be absent from various tissues of GPRC6A 2/2 mice by RT-PCR (Fig. 1C) and Western Blot (Fig. 1D). GPRC6A 2/2 mice (as well as heterozygous GPRC6A +/2 mice) were similar in gross appearance, body weight and body length to wild-type littermates (data not shown). There were no identified abnormalities in gait or physical activity between wildtype and GPRC6A 2/2 mice.

GPRC6A Deficiency is Associated with Testicular Feminization and Reduced Testosterone Levels
On closer inspection, we noted that male GPRC6A 2/2 mice had feminization of the external genitals (Figs. 2A-C). In 16-week-old male mice, the genito-anal distance (Figs. 2A and B) as well as testicular size (Fig. 2C), weight (Fig. 2D) and the weight of seminal vesicle (Fig. 2E) were significantly reduced in GPRC6A 2/2 compared to wild-type littermates. No histological abnormality of the testes was noted in GPRC6A 2/2 mice (Fig. 2F). GPRC6A was highly expressed in Leydig cells, and was also expressed in sertoli cells, spermatogonia and spermatids by in-situ hybridization analysis (Fig. 2F, lower panel). We also found abnormalities of mammary glands in male GPRC6A 2/2 mice, as evidence by greater ductal outgrowth in the mammary fat pad (in 10/14 GPRC6A 2/2 compared to 3/13 mice wild-type male mice) (Fig. 2G) and increased the mammary fat pad mass (Fig. 2H). We found no evidence of embryonic lethality or reduced fertility in homozygous null male or female mice when bred to their respective wild-type mates; however we observed a reduced litter size from breeding pairs consisting of both male and female homozygous null mice.
Next, we compared the serum testosterone and estradiol concentration in 16-week-old mice wild-type and GPR6CA 2/2 mice. Testosterone concentrations in male GPRC6A knockout mice were significantly lower and the estradiol concentrations were significantly higher in male GPRC6A 2/2 mice compared to wildtype littermates (Table 1). Estradiol levels were not different between wild-type and GPRC6A 2/2 female mice, although the circulating testosterone levels were lower in female GPRC6A null mice (Data not shown).
Since inactivation of the androgen receptor is reported to lower testosterone levels in mice [21,22,23], we examined if loss of GPRC6A lowered androgen receptor expression. We observed no difference in the amount of AR transcripts in the testis and bone marrow by RT-PCR (Fig. 3A). Real-time RT-PCR also failed to demonstrate any difference in AR expression between GPRC6A null and wild-type mice (Fig. 3B). To explore the possibility that the reduction in testosterone was due to increased aromatase- Mouse genotyping by PCR shows the presence of exon 2 in wild-type (GPRC6A +/+ ) and its absence in homozygous GPRC6A knockout (GPRC6A 2/2 ) mice. RT-PCR analysis (C) and Western Blot analysis (D) of GPRC6A expression in kidney from GPRC6A +/+ and GPRC6A 2/2 mice. RT-PCR condition, primer set and anti-mGPRC6A antibody are described in Materials and Methods. doi:10.1371/journal.pone.0003858.g001 mediated conversion for testosterone to estrogen, we examined aromatase gene (Cyp19a1) expression by real-time RT-PCR. The quantitative real-time PCR using primers located in exon 4 and the interface between exons 4 and 5 failed to document any difference in the expression of Cyp19a1 in GPRC6A null mice testis (Fig. 3C). We did, however, detect a small increase in aromatase (CPY19) protein levels in testis of GPRC6A 2/2 mice by Western Blot analysis (Fig. 3D) that was localized by immunhistochemistry ( Fig. 3E) to Leydig cells, sertoli cells, and spermatocytes site where CPY19 is known to be expressed [24]. Aromatase protein level was slightly increased by approximately 11% in testis of GPRC6A 2/2 male mice (Fig. 3D). The significance of these changes in explaining the alterations in the testosterone and estradiol levels are uncertain, since we did not measure aromatase enzyme activity.
We found no reductions in the expression of P450 17ahydroxylase gene (Cyp17), an enzyme that converts pregnenolone to dehydroandrosterone (DHEA) in the testosterone synthesis pathway [25], or estrogen sulfotransferase gene (EST/Sult1e1), which catalyzes the sulfoconjugation and inactivation of estrogens [26], in the GPRC6A 2/2 male mice testes by real-time PCR (Figs. 3G and H). Similarly, aromatase and Sult1e1 were not different in brain, fat, liver and pituitary of GPRC6A deficient and wild-type mice (data not shown). Furthermore, serum follicle-stimulating hormone (FSH) and luteinizing hormone (LH) levels in male mice were not significantly different between wild-type and the GPRC6A 2/2 male mice (Table 1). Gonadotrophin-releasing hormone gene (GnRH) message expression in the brain, however, was not altered in GPRC6A null mice (Figs. 3H and I).
Abnormalities in Renal Function in GPRC6A 2/2 mice We previously demonstrated that GPRC6A is highly expressed in the kidney [3]. We extended these observations by showing that GPRC6A is expressed in both proximal and distal tubules by in situ hybridization (Fig. 4A). Interestingly, the expression of sodiumphosphate cotransporter, NaPi IIa, (both the protein (Fig. 4B) and  16 week-old male GPRC6A +/+ and GPRC6A 2/2 mice. Data represent the mean6SEM from 6 to 10 mice in each group. * indicates a significant difference from GPRC6A +/+ and GPRC6A 2/2 mice at p,0.05, respectively. doi:10.1371/journal.pone.0003858.g002 transcript (Figs. 4C and D) was decreased ,50% in GPRC6A 2/2 mice, suggesting adaptive responses in the kidney to excrete phosphate. We also found that GPRC6A 2/2 mice had mild but significant increases in urinary calcium and phosphate excretion (calcium/creatinine ratio: 0.1960.02; phosphorus/creatinine ratio: 5.3260.31) compared to wild-type controls (calcium/creatinine ratio: 0.1360.01; phosphorus/creatinine ratio: 3.9360.28) ( Table 1). The mild hypercalciuria was not evident at 6-weeks-ofage, but was present at subsequent ages, whereas the increased urinary phosphate levels were observed only in 16-week old GPRC6A 2/2 mice. The level of serum phosphorus was also significantly higher in 16 week-old knockout mice (6.5260.18 mg/dl) compared to wild-type littermates (5.1860.21 mg/dl) ( Table 1). Circulating concentrations of calcium, PTH, FGF23, and 1,25(OH) 2 vitamin D levels were not significantly different between wild-type and GPRC6A 2/2 mice (Table 1). In addition, we found that urinary protein excretion was elevated in GPRC6A 2/2 mice by Western blot analysis (Fig. 4E). Immunohistochemical analysis revealed that this protein band in the urine represents b2-microglobulin (Fig. 4F), providing evidence for abnormalities in the proximal tubule function in GPRC6A null mice.

Metabolic Abnormalities in GPRC6A Deficient Mice
We also found that the liver of GPRC6A 2/2 mice exhibited histological features of hepatic steatosis by H&E and Oil Red O staining (Fig. 5A). Lipid positive droplets were present in hepatocytes of GPRC6A 2/2 mice but not wild-type mice. This correlated with increased triglyceride content in the livers of GPRC6A 2/2 mice (Fig 5B). In addition, we found that fasting serum glucose levels were significantly greater and insulin levels were lower in GPRC6A 2/2 mice compared to wild-type littermates (Table 1). Since fatty liver disease is a manifestation of ''metabolic syndrome'', we examined GPRC6A 2/2 mice for evidence of glucose intolerance. We performed glucose tolerance tests following IP injection of glucose (2 g/kg of body weight) after an overnight fast (GTT), and insulin tolerance tests by IP injection of insulin (0.2 units/kg of body weight) after 6 hours fast (ITT). These tests revealed that GPRC6A 2/2 mice had a significantly higher serum glucose levels during the GTT and lower sensitivity to insulin than wild-type mice in the ITT ( Fig. 5C and D, respectively).

Musculoskeletal Phenotype of GPRC6A Deficient Mice
Both male and female GPRC6A 2/2 mice had a significant reduction in lean body mass compared to wild-type littermates (7.9% and 10% in male and 11.2% and 13% in female GPRC6A 2/2 mice at 12 and 16 weeks, respectively) as assessed by PIXImus TM densitometry (Fig. 6A). There were no apparent differences, however, in muscle histology between wild-type and GPRC6A 2/2 mice (data not shown). Body fat as assessed by PIXImus TM densitometry (Fig. 6B) and white fat by gross inspection of various organs, such as testis, were increased in GPRC6A 2/2 compared to wild-type mice. Both male and female GPRC6A 2/2 mice did not have the expected age-dependent increase in bone mineral density (BMD). Indeed, BMD was significantly less at 8, 12, and 16 weeks-of-age in GPRC6A 2/2 mice as compared to age-matched wild type mice (p,0.05) (Fig. 6C). To determine if the decreased bone density might be due to defective mineralization of bone, we performed backscatter EM and bone histological analysis. Backscatter EM of cortical bone demonstrated diminished mineralization surrounding osteocyte lacunae and on the perisoteal and endosteal surfaces (Fig. 6D). Qualitative analysis of bone histology also revealed an increase in unmineralized osteoid surfaces (Fig. 6E) Serum calcium were measured at 8 week-old mice. Data are mean6SEM from more than 10 individual mice in each group. b Glucose was measured total blood after over night fasting from more than four mice of GPRC6A +/+ and GPRC6A 2/2 male mice, respectively. c Serum osterocalcin, testosterone, estradiol, LH and FSH were measured from more than six mice of GPRC6A +/+ and GPRC6A 2/2 male mice, respectively. * and ** Significant difference from GPRC6A +/+ and GPRC6A 2/2 mice at p , 0.05 and p , 0.01 respectively. doi:10.1371/journal.pone.0003858.t001 and diffuse calcein labeling of bone compared to the distinct double labels in wild-type mice (Fig. 6F), indicative of impaired mineralization.

Discussion
Recent data indicate that energy, fat, bone, and glucose metabolism are interrelated through novel hormonal regulatory pathways [20]. We observed a complex phenotype in GPRC6A 2/2 mice consisting of defective mineralization of bone and impaired osteoblast function in male and female mice, a decrease in lean body mass, an increase in fat mass, hyperphosphatemia and hypercalciuria, hyperglycemia and feminization of male mice associated with altered ratio of estradiol and testosterone, suggesting that GPRC6A participates in hormonal control of energy metabolism. Multiple organs and metabolic functions were affected by ablation of GPRC6A in mice.
Obesity-dependent metabolic syndrome is associated with impaired glucose tolerance and hepatic steatosis [27]. Similarly, GPRC6A 2/2 mice had elevated serum glucose levels (Table 1) glucose intolerance, insulin resistance and hepatic steatosis, as evidence by histological presence of fat and biochemical evidence of increase triglyceride content in livers of knockout mice (Fig. 5). The presence of a metabolic-like syndrome in GPRC6A 2/2 mice suggests that this receptor regulates metabolic pathways involved in glucose and fat metabolism, although these studies do not define the exact target organ of GPCR6A effects or whether this is a direct effect of GPRC6A or an indirect effect related to insulin resistance in GPRC6A 2/2 mice. Recently, osteocalcin produced by bone has been shown to be a circulating hormone that targets the pancreas to increase insulin secretion and adipocytes to increase adiponectin production [28]. Osteocalcin (Oc) has also been shown to target an unknown Gai coupled GPCR [29] and can activate GPRC6A in the presence of extracellular calcium [3]. These observations raise the possibility that GPRC6A, might mediate some of the actions of Oc in vivo and provide a molecular pathway linking bone and energy metabolism. Future studies will be needed to determine if this receptor mediates the effects of osteocalcin on energy metabolism and to clarify the mechanisms underlying insulin resistance and hepatic steatosis in GPRC6A 2/2 mice.
Another phenotype was feminization of male GPRC6A 2/2 mice. Low testosterone levels acting observed in GPRC6A null mice, through classical genomic mechanisms, could explain the feminization of GPRC6A null mice. The phenotype of GPRC6A null mouse resembles, but is less severe than, the AR null and Tfm (Testicular feminized) mouse models [21,23]. GPRC6A null mice, however, have inappropriately normal LH and FSH levels and increased estrogen, whereas the disruption of the genomic actions of AR is associated with increased LH and FSH and low estradiol levels [30]. Our studies do not define the mechanism of the altered testosterone/estrogen ratio. The alterations in circulating sex steroid levels in GPRC6A null mice might result from direct actions of the receptor to regulate testosterone biosynthesis or conversion to estrogens in the testis or other tissue, and/or to an indirect effect through potential GPRC6A regulation of the hypothalamus or pituitary gland. The observation that GPRC6A is expressed in the hypothalamus (GEO accession GDS565), gonadotrope cells of the anterior pituitary (GEO accession GDS2167), sertoli cells (GEO accession GDS222) and testes (GEO accession GDS2098) as well as the slight increase in level of aromatase protein expression in testis support both possibilities [2,3,14]. Testosterone replacement in GPRC6A 2/2 mice, measurement of aromatase activity and steroid biosynthesis, and assessment of GPRC6A function in other tissues will all be needed to establish the biological significance and mechanism of the observed reduction in testosterone levels in these mice.
The predominant effect of GPRC6A deficiency in bone was to impair bone mineralization. It is not clear whether these bone abnormalities represent a direct effect of GPRC6A loss from osteoblasts or secondary effects due to potential actions of the concomitant sex steroid hormone abnormalities on bone remodeling. But testosterone deficiency typically leads to increased osteoclast-mediated bone resorption [31], which was not observed in these animals. Further studies will be needed to determine the relative contribution of secondary alterations in sex hormones and primary loss of calcium-sensing receptor responses to the observed bone phenotype in GPRC6A null mice. Regardless, the expression of GPRC6A in bone and osteoblasts and the resulting bone phenotype raises the possibility that GPRC6A is a candidate for the novel osteoblastic calcium-sensing receptor [3,15]. that is distinct from CASR [15].
With regard to the kidney, we observed that GPRC6A is expressed in both proximal and distal tubules and the ablation of this receptor resulted in increased renal excretion of calcium, phosphate, and b2-microglobulin. The effects of GPRC6A loss of function on urinary calcium excretion is opposite to the hypocalcuria caused by inactivating mutations of the related calcium sensing receptor CASR [32]. It is not clear from our studies if these renal abnormalities associated with disruption of GPRC6A function are direct or indirect. We failed to observe in GRPC6A null mice evidence of secondary increase PTH or increased bone turnover that would be expected if GPRC6A ablation resulted in a renal calcium leak. Theoretically, a primary decrease in bone formation and decreased buffering capacity for calcium, with consequent increased urinary expression of dietary calcium could account for the relationship between osteopenia and hypercalciuria in GPRC6A 2/2 mice. Alternatively, a primary effect of GPRC6A to regulate gastrointestinal phosphate transport would be another possible explanation for increased serum and urinary phosphate. There are some enigmatic clinical disorders with features similar to the GPRC6A 2/2 mice that support the possibility of coordinated effects between bone formation and renal conservation of calcium. In this regard, a subset of male patients with idiopathic osteoporosis and nephrolithiasis have the combined features of decreased osteoblast-mediated bone formation and hypercalciuria, without evidence of hypogonadism, secondary hyperparathyroidism, or abnormal vitamin D levels [33] [34]. It will be interesting to determine if gene polymorphisms of GPRC6A are associated with osteopenic and hypercalciuric clinical disorders. The observed reduction of NaPi IIa in GPRC6A 2/2 mice and the low molecular weight proteinuria is consistent with a primary proximal tubular defect. The increase in serum phosphate in GPRC6A 2/2 is unexplained.
In summary, we have shown that GPRC6A has multiple functions as evidenced by abnormalities in GPRC6A null mice that include alterations in circulating testosterone and estrogen levels and feminization of male mice, defects of bone density and bone cell function and abnormalities in the renal handling of calcium and phosphate, hyperglycemia and liver steatosis. The ligand profile of GPRC6A, which includes extracellular calcium, calcimimetics, amino acids, and osteocalcin [2,3,14,20], along with the complex phenotype of GPRC6A null mice suggests that GPRC6A may represent an anabolic receptor that responds to a variety of nutritional and hormonal signals and may serve to coordinate the functions of multiple organs in response to changes of these ligands.

Generation of GPRC6A Knockout Mice
The GPRC6A-deficient mouse model was created by replacing exon 2 of the GPRC6A gene with the hygromycin resistance gene (Fig. 1). To generate the targeting construct, the hygromycin resistance gene under the control of the PGK promoter was cloned into the Sma I and Eco RV sites of pBS-lox, which was produced by cloning the oligonucleotide Lox71: 59-ctagataccgttcgtatagcatacattatacgaagttatg-39 into the Xba I and Bam HI sites and the oligonucleotide Lox 66: 59-agcttataacttcgtatagcatacattatacgaacggtag-39 into the Hind III and Sal I sites of pBluescript (Stratagene), to produce pBS-lox-PGK-Hyg. A genomic fragment from intron 1 of GPRC6A gene, representing the 59 homologous targeting region, was amplified by PCR using Advantage 2 Taq polymerase (BD Biosciences) and primers PP 232-Avr: 59-aaacctagggccattcatgaaaaaatgttgtcctcagatgaccatcc-39, and PP 233-Avr: 59-aaacctaggctcactcaacccccatgtccttccaactctagctg-39, digested with Avr II, and cloned into the Xba I site of pBS-lox-PGK-Hyg to produce pBS-GPRC6A-Intron 1. A genomic fragment from intron 2 of GPRC6A, representing the 39 homologous targeting region, was then amplified by PCR using Advantage 2 Taq polymerase and primers 298: 59-aaagtcgacctacattggtccatcgattacattagttcttgg-39 and PP 299: 59-aaagtcgacgaggccttgaggtcaaactccagaaccccagag-39, digested with Sal I, and cloned into the Sal I site of pBS-GPRC6A-Intron I to produce pGPRC6A-KO. The methods for knocking out the GPRC6A gene have been described previously in detail [35]. Briefly, the mouse embryonic cell line RW-4, derived from 129X1/SvJ mouse strain [36] and kindly provided by Stephan Teglund at the Karolinska Institute Center for Transgene Technologies, was transfected by electroporation with Not I linearized pGPRC6A-KO. Hygromycin resistant embryonic stem cell colonies were picked, expanded, and tested by PCR to identify clones in which the PGK-Hygromycin gene had correctly replaced exon 2 of the GPRC6A gene. The correctly mutated embryonic stem clones were injected into blastocysts derived from C57BL/6 3.5 days after mating and implanted into B6CBAF1 pseudopregnant females. The resulting male chimeras were bred with female C57BL/6 mice. Homozygous founders were generated by mating the resulting heterozygous mice. The successful targeting of GPRC6A in embryonic stem (ES) cells was confirmed by Southern blot analysis of the genomic DNA from ES cell clones. We observed no apparent differences in the founders generated from different ES cell clones. We focused our studies on founder line 17.
Mice were maintained and used in accordance with recommendations as described (National Research Council. 1985) Guide for the Care and Use of Laboratory Animals DHHS Publication NIH 86-23 , Institute on Laboratory Animal Resources, Rockville, MD) and following guidelines established by the University of Kansas Medical Center Institutional Animal Care and Use Committee.

PCR primers
The following intron-spanning primer sets were used for RT-PCR: mGPRC6A. 24

RT-PCR and Real-Time RT-PCR
RT-PCR was performed using two-step RNA PCR (Perkin-Elmer) as previously described [3]. Specific intron-spanning primer sets were to amply the specified transcripts. For quantitative real-time RT-PCR assessment of bone markers expression, we isolated and reverse transcribed 2.0 mg total RNA from long bone of 8-week-old mice as previously described [37]. For quantitative real time RT-PCR assessment of aromatase, CYP17 and Sutl1e1 genes expression we isolated and reverse transcribed total RNA isolated from testis, brain, fat, liver and pituitary of 33 week-old GPRC6A +/+ (n = 5) and GPRC6A 2/2 mice (n = 5) as previously described [38] using specific primer sets. For expression of AR, GnRH and NaPi IIa, the oligo primers sequences were obtained from Primer Bank, a public resource for PCR primers (pga.mgh.harvard.edu/primerbank/) [39], and the threshold cycle (Ct) of tested-gene product from the indicated genotype was normalized to the Ct for cyclophilin A as we have previously described [37].

Cell Culture and Western Blotting
HEK-293 cells were co-transfected with pcDNA3.mGPRC6A or pcDNA3 plasmid as previously described [3]. An anti-peptide antibody was raised in a rabbit against a peptide (AI-HEKMLSSDDHPRRPQIQKC) corresponding to a sequence in the extracellular domain of mouse GPRC6A (in exon 1 of mouse GPRC6A gene) produced by Abgent (San Diego, CA). For phospho-ERK, the phospho-ERK1/2 levels were determined by immunoblotting using anti-phospho-ERK1/2 mitogen-activated protein kinase antibody (Cell Signaling Technology) [3]. Urinary b2mcroglobulin was detected by rabbit polyclonal anti-b2-mcroglobulin antibody (Abcam Inc.). For aromatase expression analysis, the rabbit polyclonal anti-aromatase antibody (Abcam Inc.) and goat anti-rabbit IgG HRP secondary antibody (Santa Cruz Biotechnology Inc.) were used. Mouse anti-Actin antibody (Santa Cruz Biotechnology Inc.) was used for control protein loading.

PIXImus TM Bone Densitometer Analysis and Bone Histology
Bone mineral density (BMD) of whole skeletons and femurs were assessed at 6, 8, 12, and 16 weeks of age using a PIXImus TM bone densitometer (Lunar Corp.) as previously described [32]. Skeletons of mice were performed as reported previously by our laboratory [32].

Backscattered Scanning Electron Microscopy
Plastic embedded bone from 8 week-old wild-type and GPRC6A knockout mice were cut and polished, and mounted on aluminum sockets with bone surface facing above, sputter-coated with gold and palladium, and examined with field emission scanning electron microscopy (Philips XL30, FEI Company) equipped with a backscatter electron imaging system.

Immunohistochemistry
Wild-type and GPRC6A mouse kidney were routinely processed and embedded in paraffin. The paraffin sections at thickness of 5 mm were prepared and collected on commercially available, positively charged glass slides (Superfrost Plus, Fisher Scientific). The sections were dried on a hot plate to increase adherence to the slides. Representative sections were de-paraffined and re-hydrated through conventional methods. The sections were digested by 10 mg/ml hyaluronidase for 20 min. Nonspecific protein binding was blocked by incubation with 10% normal goat serum. The sections were incubated in polyclonal rabbit against mouse NaPi IIa (1:500 dilution) or polyclonal goat anti-human aromatase antibody (1:200 dilution) (CYP19, Santa Cruz Biotechnology, Inc.) at 4uC overnight. The negative control sections were incubated with 0.01 M PBS. Thereafter, the sections were treated sequentially with FITC-conjugated Donkey anti Rabbit IgG secondary antibody (Jackson Labs). The nucleus was stained with ready to use Hoechst (Sigma).

Whole-Mounts and Carmine Stain of Mammary Glands
Mice mammary fat pads were excised and fixed for a minimum of 2 h in Carnoy's solution (60% ethanol, 30% chloroform, and 10% glacial acetic acid). The fixed glands were washed in 70% ethanol for 15 min and then rinsed in water for 5 min. The mammary glands were stained overnight at 4uC in carmine alum stain (1 g carmine and 2.5 g aluminum potassium sulfate in 500 ml water).

In Situ Hybridization
For GPRC6A gene expression, the probe was amplified by RT-PCR using following intron spanning primer: mGPRC6A.189F (in Exon I) cgggat ccagacgaccacaaatccag and mGPRC6A.539R (spanned over Exon II and III) ccaagctt gattcataactcacctgtggc. After RT-PCR, the product was subclone into pre-cutted by BamH I and Hind III of pBluescript SK(+). Using T7 promoter will create sense RNA, and T3 promoter will create Anti-sense RNA ribo-probe.

Serum and Urine Biochemical Measurements
Serum was collected using a retroorbital bleeding technique. For urine samples collection, mice were placed in metabolic cages (Hatteras Instrument), and urine was collected for 24 h. The urine volume was measured before storage at 270uC.
Serum testosterone and estradiol levels were measured by testosterone enzyme immunoassay test kit and estradiol (E2) enzyme immunoassay test kit from BioCheck, Inc. Follicle stimulating hormone (FSH) and luteinizing hormone (LH) were measured by mouse FSH radioimmunoassay and the mouse LH sandwich assay as described by the University of Virginia Center for Research in Reproduction Ligand and Analysis Core (NICHD (SCCPRR) Grant U54-HD28934). Serum and urinary Calcium was measured by the colorimetric cresolphthalein binding method, and phosphorus was measured by the phosphomolybdate-ascorbic acid method [32]. Serum TRAP was assayed with the ELISAbased SBA Sciences mouseTRAP TM assay. Serum PTH and 1, 25 (OH) 2 vitamin D were measured the kits from Immutopics, Inc. and Immunodiagnostic system, Ltd., respectively. Serum Fgf23 levels were measured by using FGF-23 ELISA kit (Kainos Laboratories Inc.) following the manufacturer's protocol. Creatinine was measured by the colorimetric alkaline picrate method (Sigma kit 555, Sigma-Aldrich). Urinary protein and Dpd were measured by Bio-Rad and Metra Biosystems, Inc., respectively.

Metabolic Studies
For glucose tolerance test (GTT) glucose (2 g/kg body weight) was injected intraperitoneally (IP) after an overnight fast, and blood glucose was monitored using blood glucose strips and the Accu-Check glucometer (Roche) at indicated times. For insulin tolerance test (ITT) mice were fasted for 6 hr, injected IP with insulin (0.2 U/ kg body weight, Lilly Research Laboratories), and blood glucose levels were measured at indicated times as described [20]. ITT data are presented as percentage of initial blood glucose concentration.

Statistics
We evaluated differences between groups by one-way analysis of variance. All values are expressed as mean6SEM. All computations were performed using the Statgraphic statistical graphics system (STSC Inc.).


text:
A salmonid EST genomic study: genes, duplications, phylogeny and microarrays

Background Salmonids are of interest because of their relatively recent genome duplication, and their extensive use in wild fisheries and aquaculture. A comprehensive gene list and a comparison of genes in some of the different species provide valuable genomic information for one of the most widely studied groups of fish. Results 298,304 expressed sequence tags (ESTs) from Atlantic salmon (69% of the total), 11,664 chinook, 10,813 sockeye, 10,051 brook trout, 10,975 grayling, 8,630 lake whitefish, and 3,624 northern pike ESTs were obtained in this study and have been deposited into the public databases. Contigs were built and putative full-length Atlantic salmon clones have been identified. A database containing ESTs, assemblies, consensus sequences, open reading frames, gene predictions and putative annotation is available. The overall similarity between Atlantic salmon ESTs and those of rainbow trout, chinook, sockeye, brook trout, grayling, lake whitefish, northern pike and rainbow smelt is 93.4, 94.2, 94.6, 94.4, 92.5, 91.7, 89.6, and 86.2% respectively. An analysis of 78 transcript sets show Salmo as a sister group to Oncorhynchus and Salvelinus within Salmoninae, and Thymallinae as a sister group to Salmoninae and Coregoninae within Salmonidae. Extensive gene duplication is consistent with a genome duplication in the common ancestor of salmonids. Using all of the available EST data, a new expanded salmonid cDNA microarray of 32,000 features was created. Cross-species hybridizations to this cDNA microarray indicate that this resource will be useful for studies of all 68 salmonid species. Conclusion An extensive collection and analysis of salmonid RNA putative transcripts indicate that Pacific salmon, Atlantic salmon and charr are 94–96% similar while the more distant whitefish, grayling, pike and smelt are 93, 92, 89 and 86% similar to salmon. The salmonid transcriptome reveals a complex history of gene duplication that is consistent with an ancestral salmonid genome duplication hypothesis. Genome resources, including a new 32 K microarray, provide valuable new tools to study salmonids.

Background
Extensive knowledge of trout and salmon is a result of their widespread use in scientific research, as an environmental sentinel species and as a food and sport fish. Perhaps more is known about the physiology, ecology, genetics, behavior and biology of salmonids than any other fish group [1]. This background provides a wealth of data from an economically important and phylogenetically distinct group of fish that can help guide, and benefit from, new genomic studies.
The common ancestor of salmonids is purported to have experienced a whole genome duplication event between 25 and 100 MYA [6,7]. Extant salmonids are considered pseudo-tetraploid as they are in the later stages of reverting to a stable diploid state. Evidence for the ancestral salmonid autotetraploid genome duplication includes: multivalent chromosome formation during male meiosis and evidence for tetrasomic segregation at some loci [6]; one of the larger euteleost genome sizes (3-4.5 pg) with double that of sister groups Esociformes (0.8-1.8 pg, pike) and Osmeriformes (0.7 pg, smelt) [8]; homeologous chromosomal segments based on recent genetic maps and comparative studies using microsatellite markers, and duplicated gene family studies such as Hox, Major Histocompatibility complex (MH), growth hormone, and nineteen allozymes [6,[9][10][11][12].
The genome duplication in salmonids is the most recent genome duplication in this lineage. There are now a number of studies and good evidence, primarily from sequenced zebrafish and pufferfish genome sequences, for tetraploidization/rediploidization early in the ray-finned fish lineage (350-400 MYA) [13][14][15][16]. Several of these studies have suggested that the ancestral fish duplication, in addition to the two ancestral vertebrate genome duplications, are part of the reason why ray-finned fishes make up nearly half of all extant vertebrates species and exhibit tremendous biodiversity affecting their morphology, ecology, behavior and evolution.
Vertebrate species diversity and body plan diversity have commonly been linked to genome duplications, although there is some debate on how well we can draw these conclusions based on the very old genome duplications commonly studied. Mechanistically, how a genome reorganizes itself to cope with duplicated chromosomes, gene dosage effects, and the role of gene duplications for evolution and adaptation are long-standing issues in biology that remain unresolved [6,[13][14][15][16][17]. The number and diversity of salmonid species, and their relatively recent genome duplication, make salmonids ideal for examining recent events that could have played such a pivotal role in generating gene diversity and species diversity found in modern vertebrates.
The genomics resources of salmonids are being rapidly expanded through a few large-scale genomics programs [18][19][20][21][22][23]. Here we identify 354,061 new ESTs from Atlantic salmon and several other salmonid and related species in order to obtain a comprehensive view of the salmonid transcriptome, identify species relationships, identify gene duplications and introduce a new 32 K microarray tool for transcriptome analysis.

Results and discussion
cDNA libraries New, directionally cloned, mixed tissue (brain, kidney and spleen), normalized cDNA libraries were constructed for Atlantic salmon (Salmo salar; European McConnell, and Canadian, Saint John River strains), chinook salmon (Oncorhynchus tschawytscha), sockeye salmon (Oncorhynchus nerka), brook trout (Salvelinus fontinalis), lake whitefish (Coregonus clupeaformis), grayling (Thymallus thymallus), and northern pike (Esox lucius). Separate normalized libraries were constructed from Salmo salar thymus, thyroid, and head kidney tissues. In addition, one full-length, mixed tissue, large insert (> 2 kb), non-normalized library was constructed to identify longer gene transcripts. cDNA clones were isolated, purified and sequenced from the 5' and 3' ends. Clone numbers and insert sizes for the different libraries and species that were done as part of this study are listed in Table 1.

Transcript analysis: sequence and assembly
To obtain a comprehensive list of genes in salmonids, we used a strategy of deep 5' and 3' EST sequencing from a few high quality libraries. This approach complements previous studies, which examined more limited EST surveys of cDNA libraries from a large number of different tissues and developmental stages [18][19][20][21][22]. For Atlantic salmon, over 30,000 clones were sequenced from each of the thymus, thyroid, and head kidney tissue libraries. From previously described normalized libraries, [19] 9,584 additional clones were sequenced from the Atlantic salmon pyloric caecum tissue library and 60,288 additional clones were sequenced from a mixed tissue library (rgb2; Table 1). The total number of clones examined from the rgb2 library was 84,176 which yielded 127,660 sequence reads or 30% of the total Atlantic salmon EST database. Even with this deep sequencing, nearly 13% of the last 637 reads were novel (< 99% over 100 bp) and the maximum redundancy for a single transcript from the rgb2 library was 58 (Table 1).
The results of the assembly of 298,304 Atlantic salmon ESTs obtained in this study along with 138,325 ESTs from previous studies [ [19,22], GenBank] are shown in Table 2. Due to the complexities of the salmonid genome duplication and because it provides a stable, conservative starting point for all subsequent analyses, our analysis began with a first stage assembly using stringent parameters (PHRAP: 0.99 repeat stringency and 100 minscore). A second stage assembly (96% repeat stringency and 300 minscore) was implemented to combine some of contigs which may be alleles, or possibly very recent gene duplications (distinguishing among alleles, minor assembly errors, miss-calls and very recent gene duplications, particularly in lower quality sequence regions is very difficult in the absence of genomic sequence data). In Atlantic salmon, 81,398 potential transcripts (2 stage assembly) were identified, of which 29,844 (37%) were similar (BLASTX, 1e-10) to annotated sequences in CDD or SwissProt protein databases. For comparison, an assembly of 246,704 ESTs from rainbow trout [ [18,21], GenBank] resulted in 51,199 transcripts, of which 19,266 (38%) had BLASTX hits. Assembled contigs are available [24].
Transcript surveys of additional salmonid species included 4,800-7,500 clones sequenced from each of chinook salmon, sockeye salmon, brook trout, grayling and lake whitefish. 11,664 sequences were obtained from chinook salmon, 10,813 sequences from sockeye salmon, 10,051 sequences from brook trout, 10, 975 sequences from grayling and 8,630 sequences from lake whitefish. Sequence, assembly and summary statistics are shown for those data obtained in this study (Table 1) and when combined with data from public databases ( Table 2). In addition, to provide non-genome-duplicated sister group comparisons, 2,304 clones were sequenced from northern pike (3,624 sequences) ( Table 1 and 2). For many of these species, the ESTs provided in this study represent nearly all or most of the known transcripts. Recently published data from rainbow smelt (Osmerus mordax) [23] was also included in Table 2.
To examine the relationships among the contig consensus sequences of Atlantic salmon we compared all contigs (including singletons) against each other by BLAST and plotted the number of top pair-wise alignments (E-value < 1e-50; length > 200 bp) with the identity score ( Figure  1). 36,775 contigs showed greater than 80% identity over 200 bp to at least one other contig. Of these, 12,883 were 97-99.9% similar to at least one other contig. These contigs may represent alleles, recent duplicates or errors in sequence data. 23,892 contigs show between 80 and 96.9% identity with at least one other contig. The large number of duplicated transcripts observed in the Atlantic salmon genome is consistent with the hypothesis of an ancestral salmonid genome duplication, though it is surprising that so many of the duplicated contigs are so similar. This observation is being pursued further in a separate study. The analysis of contig similarity shows that the majority of the 81,398 contigs represent distinct transcripts. Note that since the assembly process itself combines sequences with high levels of similarity (> 96% repeat stringency with minscore > 300; see Methods), very recent duplications may not all be identified in this process. Furthermore, since the species used in this study differ by greater than 5% (Table 3), this process would be expected to identify ancestral salmonid duplications occurring at or prior to the rainbow trout and Atlantic salmon speciation.
Determining the number of genes in Atlantic salmon from the number of EST contigs is difficult for several reasons; 1) the partial representation of genes by EST sequences may result in several contigs associated with a single gene transcript, 2) allelic or recently duplicated genes may be represented by similar but unique transcripts (this latter case is particularly problematic in pseudotetraploid salmonids), 3) alternative splicing, alternative poly adenylation and termination sites from the same gene can result in different transcripts, and 4) transcription products can occur from intergenic regions. An estimation of the number of genes in salmonids will require additional information such as full-length cDNA sequences and gene mapping information.

Salmonid comparisons
Similarity among the different salmonid species was assessed using the top BLASTN hit against Atlantic salmon and rainbow trout EST contig databases. The similarity values from chinook salmon, sockeye salmon, rainbow trout, Atlantic salmon (McConnell and Saint John River strains), brook trout, grayling, lake whitefish, northern pike and rainbow smelt are shown in Table 3   Number of aligned contigs (y-axis) out of 81,398 total contigs is plotted against percent similarity of alignments (x-axis) Figure 1 Number of aligned contigs (y-axis) out of 81,398 total contigs is plotted against percent similarity of alignments (x-axis).
contigs matching (E-value < 1e-25) at least one contig in the rainbow trout or Atlantic salmon databases. These comparisons provide only a very general indication of the similarity between transcriptomes of various salmonids, as assemblies contain both 5' (generally genic regions) and 3' (generally 3'-UTR regions) transcript reads. However, these DNA sequence similarity values correspond well to the limited number of values in the literature. Non-coding sequence similarity between rainbow trout and Atlantic salmon are 95% over 120 kb in MH class IA and B loci [12], and 93-97% over 4 kb in growth hormone (GH) genes [11]. Similarity between salmon and whitefish is 90-93% in GH genes [11].
Northern pike and rainbow smelt average 89.4 and 86.1% identity to rainbow trout and 89.6 and 86.2% identity to Atlantic salmon, but only 25-39% of these contigs matched anything in the rainbow trout or Atlantic salmon database. These latter comparisons have many fewer significant similarities identified partly because of the much older divergence times [3]. However, the reason for the lower than expected number of matches between northern pike and rainbow trout or Atlantic salmon is not clear. While the more distantly related rainbow smelt contigs show similar numbers of BLASTX hits to protein databases as salmonids, the northern pike contigs showed very few similarities to Atlantic salmon and rainbow trout contigs (25% compared to 39% for rainbow smelt and 87% for lake whitefish) and very few BLASTX hits to protein databases (8% compared to 50% for rainbow smelt and 36% for lake whitefish). One possible explanation may be due to longer 3'-UTRs in northern pike, but this remains to be confirmed.

Transcriptome representation
It is difficult to assess how comprehensive the extensive Atlantic salmon and rainbow trout EST databases are. However, 73% (37,573 of the 51,199) of all rainbow trout contigs are also found in Atlantic salmon. Moreover, only 28% of those transcripts unique to rainbow trout (13,626) have protein hits (E < 1e-25) that support their legitimacy as genic regions, while other single ESTs may be from spurious transcription. 91% of lake whitefish transcripts have a significant similarity (BLASTN comparisons with e-values less than 1e-25) to the Atlantic salmon or rainbow trout databases. Comparative data from chinook salmon, sockeye salmon, brook trout, grayling, lake whitefish and rainbow smelt are provided in Table 3.
Overall, these data provide support for extensive gene coverage in salmonid EST databases.

Full-length analysis
The rapid progress of EST sequencing has enabled an estimation of the number of full-length cDNA clones. Fulllength cDNAs (fl-cDNAs) are defined as having a "Start -Open Reading Frame (ORF) -Stop -3' UTR -polyA signal" with the ORF corresponding to a full-length protein.
Given multiple start and stop sites, alternative splicing and partial homologies to known proteins, it is difficult to give precise numbers of completed fl-cDNAs. However, TargetIdentifier (using BLAST comparisons to full-length genes in databases and Start signals; [25]) identifies 17,399 possible fl-cDNAs (averaging 1,361 bp in length) from the 81,398 possible transcripts in Atlantic salmon and 10,453 fl-cDNAs from the 51,199 rainbow trout transcripts. Thus far, about half of the predicted fl-cDNA meet all of the criteria above, and many of the fl-cDNAs are already fully characterized on a single clone. These tend to be the shorter (< 1.5 kb) genes. The list of over 10,000 putative fl-cDNA transcripts assembled from ESTs is available at the GRASP website [24] and further identification of clones for complete sequence analysis is underway.

Salmonid EST, assembly, ORF and annotation database
All ESTs have been deposited in GenBank, however the EST assemblies themselves and the resulting consensus sequences are also very useful in identifying genes. These assemblies, together with the raw data are available [24]. The assembly consensus sequences are available for download and for searching using BLAST tools. A contig visualization tool was developed to allow users to search for similar consensus sequences using BLAST searches, identifying consensus names and then visualizing the sequences, alignment, open-reading frames (ORFs), Tar-getIdentifier predictions, and BLASTX hits in a single view ( Figure 2: Cluster tools). Until such time as the genomes are completed, this database provides the salmonid community with access to several levels of EST and gene analyses.

Salmonid phylogeny and gene duplication
The relationships among major groups of salmonids have been largely unresolved, particularly with respect to the placement of Salvelinus (represented in this study by brook trout), Oncorhynchus (represented here by rainbow trout, chinook and coho salmon) and Salmo (Atlantic salmon) within Salmoninae, and the placement of Thymallinae (grayling), Coregoninae (whitefish) and Salmoninae (salmon) within Salmonidae [2,11,[26][27][28][29][30]. From the EST contigs (Table 2), 78 separate gene sets have been identified, each of which contained at least one EST contig sequence from each of five major salmonid genera (Oncorhynchus, Salmo, Salvelinus, Coregonus, and Thymallus), in addition to representation by a non-salmonid (Osmerus). Contig sequences within each gene set were aligned, trimmed to a common length (minimum of 300 bp) and analyzed using phylogenetic methods. 73 of the 78 gene sets could be identified by BLASTX searches to SwissProt databases (Table 4). For each gene set, a 70% neighbour-joining (NJ) consensus tree based on 500 Screen shot of Atlantic salmon contig viewer bootstrap replicates was generated and the consensus tree rooted with Osmerus mordax sequences (rainbow smelt).
The single species tree shown in Figure 3 represents a compilation of the phylogenetic results from 78 gene sets. In the summary tree, each branch is noted by; i) the number of 70% consensus NJ trees supporting the branch ii) the number of 70% consensus trees providing no resolution to the branch point, and iii) the number of consensus trees that conflict with the shown result. In this summary, the placement of Salmo as a sister group to Oncorhynchus and Salvelinus is supported in 18 of the 27 gene consensus trees for which resolution was found. Eight alternative consensus trees support grouping Salmo and Salvelinus, one consensus tree supports grouping Salmo and Oncorhychus, and the remaining 51 trees provide no resolution. Thus the overall result is in agreement with some of the more recent studies examining mitochondrial and nine nuclear genes [26], and suggests good support for grouping Oncorhynchus and Salvelinus apart from Salmo within the Salmoninae subfamily.
Consistent with traditional nomenclature, the Salmoninae group, which includes Salvelinus, Oncorhynchus and Salmo is also very well supported with 51 of the 54 resolved trees consistent with this grouping. The three discrepant trees supported a Salmo/Coregonus grouping.
The relationships among the three subfamilies within Salmonidae have not been extensively addressed at the molecular level. However, on the basis of a morphological analysis, Coregoninae (whitefish and ciscos) has been hypothesized as the earliest branch within the salmonids [ [29], also see [30]]. In the present analysis, 14 of the 25 informative gene sets are more consistent with the basal position of Thymallinae ( Figure 3). Of the discrepant trees, 8 sets support a Thymallinae/Coregoninae grouping and 3 support an ancestral position of Coregoninae. While these data are not definitive, there appears to be some support for an ancestral Thymallinae branching within the Salmonidae with Coregoninae as the sister group to Salmoninae. These data provide the first largescale molecular view of salmonid subfamily relationships and provide an important perspective on future analyses of duplicated genes, as well as physiological and ecological traits [27] that have evolved subsequent to the ancestral salmonid genome duplication.
The salmonid whole genome duplication hypothesis makes it is difficult to separate an analysis of species relationships from gene phylogeny. One expectation arising from a relatively recent genome duplication is evidence for extensive nuclear gene duplicates. Subsequent to the genome duplication, the number of observed duplicated transcribed genes is expected to decrease as, over time, one of the duplicates becomes transcriptionally inactive.
When multiple species are examined, some species may have both duplication products while other species may have only one representative. Evidence of an ancestral duplication is identified in gene trees that contain multiple species trees that may have missing representatives. Of the 78 gene sets examined in this study, 51 show clear evidence of multiple species trees within gene trees that are consistent with a gene duplication in the ancestor of Salmonidae, sometime after the separation of Osmeriformes and Salmoniformes fish. 23 gene sets (Table 4) provided no evidence for any ancestral gene duplication, and 4 sets could not be interpreted. The data from 78 gene sets representing 372 consensus sequences and 11,397 bp of aligned DNA from five salmonid genera, indicate that a large number of salmonid genes show evidence of extensive gene duplication at a phylogenetic position that is consistent with the whole genome duplication in the ancestral Salmonidae hypothesis. Further studies of Esociformes fish will more precisely establish the timing of some of these gene duplications.

Salmonid 32 K microarray
To use the data generated by ESTs and assemblies for examining gene expression, a new 32,000 feature cDNA microarray was developed. This new array is based on the existing 16 K GRASP array [31] plus 14,496 additional Atlantic salmon and 1,491 additional rainbow trout contigs that were identified as unique and were successfully amplified in this study. The 32 K cDNA microarray is composed mainly of 27,917 Atlantic salmon (AS) and 4,065 rainbow trout (RT) cDNA elements or features. 54% of the elements have fairly stringent (1-e10) hits to annotated members in public protein databases. Hybridization performance of this array was evaluated using Atlantic salmon, rainbow trout, coho salmon, brook trout and lake whitefish RNA obtained from liver organs. The success of hybridization of labeled target to the salmonid elements was judged by the numbers of Atlantic salmon and rainbow trout elements passing background plus 2 SD threshold values (see Methods). No transformations or normalizations were performed on the data. Overall statistics are presented in Table 5. In summary, for RNA isolated from the liver of Atlantic salmon, rainbow trout, coho salmon, brook trout and lake whitefish, an average of 48% of the 32,018 elements showed significant detection levels of expression. Comparing these results to that from the previous 16 K GRASP arrays indicates that doubling the number of elements from 16 K to 32 K resulted in the ability to assess expression patterns of approximately 61% additional transcripts. This represents a substantial increase in our ability to assess gene transcription patterns in salmonids. The hybridization performances of the different salmonid species (assessed from numbers of Atlantic salmon and rainbow trout elements passing threshold) conformed to expectations, given the close  evolutionary relationships of the species tested (92-94% identity, Fig. 3) and, with the possible exception of brook trout, all members of the family Salmonidae tested showed similar levels of hybridization to the Atlantic salmon and rainbow trout elements on the 32 K microarray. As the Salmonidae family represents 68 closely related species, the 32 K cDNA array provides an excellent opportunity to evaluate gene expression patterns of a large group of culturally and economically important species.

Conclusion
Atlantic salmon and rainbow trout now rank 19 th and 29 th in terms of species representation in EST databases with over 730,000 salmonid ESTs in total. Almost half of these data are presented in this study. These data provide an excellent genetic resource for physiological, ecological, biochemical, behavioral, disease and biological studies of salmonids. They also provide key materials for the development of polymorphic markers for genetic and physical genomics maps, for the identification and analysis of proteins and for the development of microarrays and primers for transcriptional analyses.
Summary of 78 gene set consensus (70%) trees depicting the relationships among the major groups of Salmonidae  Table 4.  * Percent elements on cDNA array with median signal intensity greater than threshold (background signal+ 2SD). %CV is percent coefficient of variation and "n" is the number of biological replicates.

Subfamily
Transcript assemblies and analyses have identified over 81,000 possible transcripts from Atlantic salmon and 51,000 transcripts from rainbow trout. These assemblies and consensus sequences are available from the author or through a database housed on the GRASP website [24]. As many as 17,399 full-length Salmo salar gene assemblies are present in this database.
Comparison of orthologous ESTs from Atlantic salmon, rainbow trout, chinook salmon, sockeye salmon, brook trout, lake whitefish, grayling, northern pike and rainbow smelt show that Pacific salmon (Oncorhynchus), Atlantic salmon (Salmo salar) and brook trout (Salvelinus fontinalus) average 94-96% similarity. Lake whitefish (Coregonus clupeaformis) and grayling (Thymallus thymallus) are more distant from Pacific and Atlantic salmon (92%), followed by northern pike (Esox; 89%) and rainbow smelt (Osmerus; 86%). A view of salmonid relationships and support for the salmonid genome duplication has been found. With the new EST database, a new, more extensive 32 K cDNA microarray has been developed to help assess gene expression patterns in salmonids.

Tissues, RNA, Aquaculture and Sampling
Salmo salar (McConnell strain), Oncorhynchus tshawytscha and Oncorhynchus nerka tissues were obtained from the Department of Fisheries and Oceans (Robert Devlin, WestVan Lab., West Vancouver, British Columbia). Salvelinus fontinalis and Coregonus clupeaformis tissues were obtained from Louis Bernatchez (Laval University, Quebec). S. salar (Saint John River strain; brain, kidney and spleen) were obtained from Vanya Ewart (NRC Institute for Marine Biosciences, Nova Scotia). Thymallus thymallus brain, kidney and spleen tissues were obtained from Craig Primmer (University of Turku, Finland). Esox lucius were captured by gill net from Charlie Lake British Columbia. All fish were euthanized, followed by rapid dissection of tissues. Tissues were flash frozen in liquid nitrogen or dry ice and stored at -80°C until RNA extraction.

cDNA libraries
Total RNA or poly(A)+ RNA (FastTrack MAG kit; Invitrogen) was extracted from flash frozen tissues. Salmo salar and Oncorhynchus tshawytshca mixed tissue (spleen, head kidney, brain) libraries were directionally constructed in both pCMV Sport-6.1 (Research Genetics Inc.) and pAL-17.3 (Evrogen Co.). S. salar (normalized head kidney, thymus and thyroid), Coregonus clupeaformis. Thymallus thymallus and Salvelinus fontinalis libraries were constructed in pAL-17.3 (Evrogen). The Oncorhynchus nerka mixed tissue normalized library was also constructed in pCMV Sport-6.1 (ResGen). S. salar (mixed tissue St. John strain) and Esox lucius libraries were constructed in pDNR-Lib using Creator SMART cDNA library construction kits (Clontech). Insert sizes of cDNA libraries were determined by visual comparison of clone restriction fragments with the DNA size markers HindIII (GibcoBRL) and 1 kb ladder (GibcoBRL).
The number of Salmo salar contigs was assessed using the PHRAP assembly program because of its ability to assemble very large numbers of ESTs in a single run, and its integration with PHRED base quality scores on primary reads and subsequent consensus sequences. The CAP3 assembler [39] was also used and similar results were obtained for smaller datasets. For this study, contig assembly employed a two-stage process. The first stage assembly used parameters 100 minscore and 0.99 repeat stringency to build contigs and consensus sequences that appeared to separate alleles of many transcripts. The second stage used the consensus sequences (with quality scores) from the first stage and parameters 96% repeat frequency and 300 minscore to build contigs and consensus sequences that appeared to combine some of the contigs that contained some base calling discrepancies, as well as what appeared to be alleles or very recently duplicated genes. Various parameters were tested and final parameters were chosen to minimize the number of contigs, where the number of contigs changed the least with respect to small changes in parameter values, and where distinct contigs appeared to have some biological significance (i.e., 99/ 100 appeared to separate many alleles and 96/300 as a second stage appeared to join some alleles and provided values that separated a clear majority of orthologous salmonid gene comparisons). With both sets of parameters, we were able to discriminate between similar sequences from different salmonid species. Sequences in contigs containing more than one polyA site were removed from the assemblies as they may represent chimeric clones.
Assemblies provide rough estimates of transcripts. Several algorithms have been examined and all have strengths and weaknesses. Examples of other assemblies include DFCI gene indexes [40,41] that estimate 83,554 TCs+singletons from 244,984 rainbow trout ESTs and 63,138 contigs from a partial 236,009 EST dataset from Atlantic salmon (these assemblies are periodically updated). INRA [21] using CAP3 estimates 56,392 transcripts (contigs + singlets) from 326,719 rainbow trout ESTs and 45,349 contigs from a partial Atlantic salmon EST database. Uni-Gene [42], from NCBI does not provide true assemblies and may cluster duplicated genes into single bins, which is problematic in salmonids. UniGene estimates approximately 30,000 and 25,000 UniGene sets in Atlantic salmon and rainbow trout respectively. While differences exist, the general number of estimated transcripts is similar. Problem areas that have been identified in assemblies tend to be associated with long transcripts, so these contigs will have to be treated carefully and perhaps manually edited. Assemblies are freely available from the author and the GRASP website [24]. As a caveat, because of the purported recent duplication of the salmonid genome and potential for miss-assembly of duplicated transcripts, these contigs have to be treated with caution.
Percent identity measures between contig consensus sequences from the various species were obtained from BLASTN alignments where a minimum length of 200 bp was observed. As in other distance measures, this finds the most similar sequence fragments and is biased high, particularly for more distant comparisons. A partial estimate of the impact on more distantly related sequence comparisons is the increased number of contigs for which no cross-species alignments were found and the reduction in average length of alignments. These values are provided in Table 3. However, the percent identity measure provides an estimate of observed similarity that is useful for evaluating potential cross-species DNA hybridizations in microarray experiments (see below).

Gene phylogenetic analysis
Contig sequences from Salmo salar (Atlantic salmon), Oncorhynchus mykiss (rainbow trout), Osmerus mordax (rainbow smelt), Coregonus clupeaformis (lake whitefish), Salvelinus fontinalis (brook trout), and Thymallus thymallus (grayling) [ Table 2] were BLASTed against each other (evalue < 1e-35, hits > 100bp) and the results used to gen-erate clusters of contigs. Bins of similar sequences, or clusters, were generated containing all contigs irrespective of species origin that had alignments with greater than 75% of the length of the shorter sequence and had greater than 70% identity in the overlapping regions (alignments consisted of ends-free alignment with scores of 2/-2/-5/-1 for match/mismatch/open gaps/extend gaps [43]. After the contigs had been grouped into clusters, the individual clusters were then further selected to only contain contigs that had mutually overlapping regions and all contig members were trimmed to the largest common alignment (same alignment parameters as above). A good alignment was considered to be greater than or equal to 300 bp in length with greater than 60% identity in the overlapping region. At this point, clusters that did not contain at least one sequence from each of the six target species were discarded. This resulted in a dataset of 78 clusters or gene sets. All gaps (and their corresponding positions in other sequences of the cluster) were removed, and the data within each gene set were bootstrapped 500 times. The PHYLIP package was used because it offers many different analysis methods, is freely available and is commonly used [44]. Distance matrices were computed for each bootstrapped dataset within each cluster using the F84 model of nucleotide substitution and Gamma-distributed rates of variation across sites with a coefficient of variation of 0.5 [44]. Neighbor-joining trees were then computed from each set of distance matrices and the set of resulting bootstrapped trees was used to derive a 70%-majority consensus tree [44]. The consensus trees were rooted with Osmerus mordax, and simplified by iteratively collapsing all pairs of leaf nodes having the same species and showing > = 98% similarity in the aligned portion of their sequences. Independently, maximum likelihood trees were generated for all 78 data sets using the default options with the Phylip program dnaml (transition/transversion ratio of 2.0, empirical base frequencies, constant rate variation among sites). A general evolutionary model was used for the 78 data sets because each set potentially consisted of a mixture of unidentified coding and noncoding data. All of the 78 ML trees were consistent with their 70%-consensus bootstrapped Neighbor-joining counterparts. EST accession numbers used to make contig consensus sequences, alignments and the 70% consensus trees are available [see Additional file 1] or online at the GRASP website [24].

Microarray Clone selection
Starting from the existing GRASP 16 K cDNA microarray [24], additional clones were selected for representation on the following basis: a) the contig ( Table 2) includes at least one clone that is on hand; b) the contig is of high quality with few conflicting positions, few singleton positions, no interior singleton positions (potential chimeric sites) and there are at least two clones in the contig (from at least 2 plates, and preferably from at least 2 libraries); c) if the contig is singleton then it must have a good BlastX hit (e-value < 1e-8) or other indication of orientation (eg. consistent poly(A) tail information); d) contig must have < = 94% identity to another sequence on the chip (the existing 16 K plus any new contig; not counting rainbow trout orthologs); and e) no tRNA, ribosomal, or mitochondrial sequences. We chose clone representatives within each contig based on: a) the reliability of the cDNA library and sequence; b) high similarity to consensus of contig (allow 20 bp at ends for poor trimming); c) reliable sequence from the 3'-end of contig and correct (3' -> 5') orientation; and d) ownership of clone.

Microarray fabrication
The initial clones were robotically rearrayed from daughter glycerol stock 384-well plates into 96-well plates prefilled with 8% glycerol in 2XYT + ampicillin with a MicroGrid II-610 (Biorobotics, Cambridge, UK), incubated overnight at 37°C, and checked for uniform optical density. Plasmid inserts were PCR-amplified in a MJ Tetrad PTC-205 thermocycler (Bio-Rad, Hercules, CA, USA) by using 1.0 μL overnight culture, 0.3 μM M13/pUC forward primer (5'-CCCAGTCACGACGTTGTAAAACG-3'), 0.3 μM M13/pUC reverse primer (5'-AGCGGA-TAACAATTTCACACAGG-3'), 2 mM MgCl2, 10 mM Tris-HCl, 50 mM KCl, 200 μM dNTPs, 1U AmpliTaq (Roche Diagnostics, NJ, USA), and nuclease-free H2O (Qiagen, Valencia, CA, USA) to 100 μL. PCR conditions were as follows: 2 min at 95°C denaturation; 35 cycles of 30 sec at 95°C, 45 sec at 59°C, and 4 min at 72°C; and 7 min at 72°C. Hotstar taq (Qiagen) was used to amplify additional inserts (clone set 2) with an initial denaturation of 15 mins. Amplicon specificity and yield was analyzed by capillary electrophoresis using the HT DNA SE 30 LabChip on Caliper AMS 90 system (Zymark-Caliper Life Sciences, MA, USA). PCR products were robotically cleaned (Qiagen) and consolidated into 384-well plates, lyophilized by speed-vac, and resuspended in 20 μL 3× SSC plus 1.0 M betaine. All cDNAs (average printing concentration of 165 ng/ul [original inserts] and 100 ng/ul [new inserts]) were printed as single spots on Erie Aminosilane slides (Erie, Portsmouth, N.H., USA) with a Genetix QArraymax microarray printer (Genetix, New Milton, Hampshire, UK) or MicroGridII-610 printer (Biorobotics, Cambridge, UK). All clones and controls were distributed randomly on the array. Genetix aQu 65 um quill pins or Biorobotics 10 k quill pins in a 48-pin tool were used to deposit < 1.0 nL (0.1 ng cDNA) per spot onto the slide. The resulting microarrays have a 4 × 12 subgrid layout with 699 spots per subgrid, each spot having diameter and pitch of 90-130 and 160-190 μm, respectively. A 280-bp GFP (green fluorescent protein) cDNA was amplified from a GFP clone (BD Biosciences, Mountain View, CA, USA) by using the primers (5'-GAAA-CATTCTTGGACACAAATTGG-3') and (5'-GCAGCTGTTA-CAAACTCAAGAAGG-3'), and printed in subgrid corners to assist in placing on the grid. The slides were crosslinked in a UV Stratalinker 2400 (Stratagene, La Jolla, CA, USA) at 300 mJ. One slide every 20 to 30 slides was hybridized with labeled random 9-mer oligonucleotide (SpotQC, Integrated DNA Technologies, Coraville, IA, USA) and scanned using GenePix 4200AL scanner (Molecular Devices, Sunnyvale, CA, USA). Presence/absence, shape, signal intensity vs. background, diameter and DNA binding site capability were measured for each spot on the slide using files generated by Imagene software (BioDiscovery Inc., El Secundo, CA, USA). Position and description of flagged spots (spots absent or thought to be unusable during post hybridization analysis), sub-grid defects and other noticed irregularities are recorded. Two PCR fragments from each plate were randomly selected and sequenced to ensure correct matches to the original clone sequence in the EST database. For controls, Stratagene SpotReport Alien cDNA Array Validation system PCR products (Cat # 252550) composed of 10 unique PCR products are spotted five times on the array. Corresponding mRNA for these PCR products can be purchased from Stratagene. The alien mRNA spikes can be used to determine mRNA quality, cDNA synthesis efficiency, positive and negative hybridization control, normalization for dye differences and determination of hybridization consistency.

Microarray hybridizations
The microarray experiments were designed to comply with MIAME guidelines. To minimize technical variability, all targets were synthesized in one round and hybridization experiments were conducted on slides from a single batch. Each hybridization experiment included dye-flips to compensate for cyanine fluor effects. Total RNA samples were quantified and quality-checked by spectrophotometer and agarose gel, respectively. All hybridization experiments were performed using the SuperScript III Indirect cDNA Labeling System kit and following manufacturers instructions (Invitrogen). Briefly, total RNA was reverse transcribed using an anchored oligo d(T) 20 primer in cDNA synthesis reactions that incorporated aminoallyl-and aminohexyl-modified nucleotides. The modified cDNAs were then labeled with fluorescent Cy5 or Cy3 dye in reactions with the amino-functional groups in coupling buffer.

Microarray analyses
Fluorescent images of hybridized arrays were acquired immediately at 10 um resolution using ScanArray Express scanner (PerkinElmer). The Cy3 and Cy5 cyanine fluors were excited at 543 nm and 633 nm, respectively, at the same laser power (90%), with adjusted photomultiplier tube settings between slides to balance the Cy5 and Cy3 channels. Fluorescent intensity data was extracted from TIFF images using Imagene 5.6.1 software (Biodiscovery). Quality statistics were compiled in Excel from raw Imagene fluorescence intensity report files. The hybridization performance of labeled targets to salmonid features was assessed as a percentage of features bound from the numbers of AS and RT features passing a hybridization signal threshold. Signal threshold was defined by 2 standard deviations above the signal mean for the 3× SSC/ betaine buffer spots. Outliers of buffer spots were removed based on the Median Absolute Deviation method [45] whereby elements with a test statistic value greater than 5 were removed. No transformations or normalizations were performed on these data. Only features deemed present by Imagene 5.6.1 (excluding marginal and absent values) were used for analyses.