Four major points emerge from the results obtained in the present study. The first one regards the high level of flexibility, reliability, and reproducibility of the microarray platform developed for the European sea bass. Oligonucleotide probe design proved to be feasible on nearly all available transcribed sequences and almost every designed probe yielded a successful hybridization signal, providing a most complete representation of the current transcriptome for D.labrax
. One of the most important requirements for a microarray platform is good reproducibility. Analysis of biological replicates showed extremely high correlation coefficients, demonstrating great reproducibility of microarray data. As already reported [16
], the use of a one color labeling scheme has apparently no negative effects on the quality of data, but it greatly simplifies the experimental design and allows easier analysis of novel samples. Finally, validation with qRT-PCR, a method that is based on a different technical approach, confirmed microarray hybridization data with the exception of a single gene, CDK2-associated protein 1, which showed a positive, but not significant correlation between results of the two different methodologies (Table ). This is likely due to the small difference in expression for this gene (mean fold change estimated from array data is 1.3). Lack of correlation between microarray and qRT-PCR for genes exhibiting low levels of change (< 1.4 fold) has been commonly reported. Indeed a two-fold change is usually considered as the cut-off below which microarray and qRT-PCR data begin to lose correlation [26
The second point concerns the annotation of the sea bass transcriptome. Whilst the design and validation of a custom oligo DNA microarray has become a straightforward procedure, a proper functional classification for the collection of unique transcripts that are represented onto the DNA microarray is certainly more difficult to achieve in non-model species, particularly in teleost fish.
High-throughput technologies have now prompted different comprehensive approaches toward functional analysis of lists of differentially expressed genes/proteins, Gene Ontology [27
] being the most popular of them. These methods offer relevant advantages, although gene ontology annotations suffer from known limitations [28
]. In non-model species these approaches are based on putative homology against better characterized organisms rather than on direct experimental evidence. There are here two major sources of error; the first related to the problem of identifying true orthology [29
], the second to the fact that sequence homology does not always imply functional homology, due to inherent differences between the species compared. In the case of teleost fish, extracting gene functional annotations from mammalian models poses an additional problem. The ancestral teleost genome underwent a whole genome duplication (WGD) after the separation from the tetrapod lineage [30
], which led to the presence of two copies of each gene present in higher vertebrates, a fraction of which genes was retained through different mechanisms (e.g.
neofunctionalization, subfunctionalization, genomic inertia), whereas others have been lost by gene deletion or pseudogenization [31
]. A ready example of such phenomenon is provided by the gene that showed the largest fold change difference in the present study, Arrestin 3 (Retinal X arrestin, ARR3) (see Additional file 3
). ARR3 belongs to the cluster of beta-arrestins, which includes four paralogs (Arrestin beta 1, ARRB1; Arrestin beta 2, ARRB2, S-AG arrestin, ARR3) in the human genome [32
] but at least seven paralogs (ARRB1a; ARRB1b, ARRB2a, ARRB2b, S-AG arrestin, ARR3a, ARR3b) in the stickleback genome. Despite all these limitations, functional annotation through sequence homology remains the best available method in non-model species. In the case of D. labrax
, a good percentage (>60%) of DLPD entries has a significant match with a known gene/protein, generally higher than the value observed in similar studies on other species (gilthead sea bream 40% [16
], turbot 50.7% [33
] Senegal sole 40.6% [34
], largemouth bass 46% [10
], pre-smolt Atlantic salmon 50.3% [35
], channel catfish 51% [36
], and comparable to that reported for the Atlantic halibut 60% [37
]). A smaller proportion (>4,500) of sea bass transcripts could be associated with a GO term, most likely as a result of the more stringent criteria enforced in Blast2GO and the lack of GO annotation for part of the protein sequences matching DLPD entries. Nevertheless, sufficient information was obtained to construct a custom background, which takes into account the fact that the sea bass transcriptome is only partly represented onto the microarray, and to identify biological processes that are significantly enriched among differentially expressed genes (see Table ). A second approach of annotation by similarity, through direct or indirect comparison against a model species transcriptome appears to convey comparable information, with over 6,000 human or 5,000 zebrafish putative homologs. The use of the human transcriptome as a reference provides a larger number of significant functional annotations (see Table ), which is somehow expected since the probability of identifying significantly enriched GO terms depends in part on the size of the gene list to be analyzed and the corresponding background. On the other hand, using human-centered annotations might be more prone to the risks described previously. Finally, the results obtained with different methods (GOStat compared to DAVID, DAVID human knowledgebase against zebrafish one) show a certain degree of overlap, yet there are terms that appear only in one or two analyses, therefore the obtained results should be interpreted with some caution.
While the majority of DLPD entries can be associated to a known gene, a substantial "silent" minority exists of transcripts that could not find any significant match against a variety of sequence databases. Even if not associated with a known coding gene, these transcripts might convey useful information. For instance, a large proportion (190/242 with a threshold of e-5, and 169/242 with a threshold of e-10) of all differentially expressed transcripts in jaw deformed sea bass find a significant match against a specific region in the stickleback genome. The strong conservation of large genome segments between D. labrax and G. aculeatus (R. Reinhardt, unpublished data) and the availability of a physical map of the sea bass genome (F. Galibert, personal communication) will allow the use of matching DLPD transcripts as positional candidates to identify loci involved in the genetic determination of mandibular prognathism.
The 40% of non-matching DLPDs likely represents different types of transcripts, e.g.
5' and 3' end untranslated regions (UTRs) or alternative splicing isoforms that have not been (yet) characterized in other species, extremely fast-evolving protein-coding regions, novel genes. In addition, a fraction of non-matching sea bass transcripts might belong to the expanding universe of non-coding RNAs, which appear to cover an increasing part of the animal genome [38
]. Although originally dismissed as transcriptional noise, evidence is accumulating for a functional role of these transcripts [40
]. The third point of the present study, the presence of NATs in the sea bass transcriptome, is related to this issue. Natural antisense transcripts have been originally identified searching EST collections, and appear to be widespread across animal species, albeit at diverse frequency [41
]. Various putative functions have been proposed for NATs [42
], with an increasingly relevant role in the production of endogenous siRNAs [43
]. Oligo DNA microarrays have been used to specifically detect sense-antisense gene expression [44
], although there are some caveats about the risk of false positive NAT detection due to genomic DNA contamination of RNA extracts or unintended labeling of both cDNA strands [46
]. Microarray analysis of the sea bass transcriptome was not specifically aimed at investigating NATs yet it provided preliminary evidence for the existence of putative antisense RNAs in this species. Based on sequence homology with the stickleback genome, in most cases sea bass NATs represent non-overlapping regions of corresponding sense transcripts (see Additional File 6
), and can derive from protein-coding regions as well as non-protein-coding sequences, including exonic, intronic and intergenic sequences. While the biological roles of NATs remain to be elucidated, putative sea bass NATs show lower levels of expression compared to the corresponding sense DLPD transcripts (see Additional File 7
), as already reported for other species [46
]. However, more than 10% of putative NATs showed higher expression levels compared to sense transcripts (e.g.
14 in 38 days-old larvae with fold change ranging between 1.4 and 71). Functional analysis of sea bass NATs also showed significant enrichment of certain molecular functions, although these results require further confirmation.
The fourth and last evidence from the present study was the discovery of a set of differentially regulated transcripts in the mandible of jaw-deformed juvenile sea bass compared to normally developed animals. Developmental defects are generally thought to trace back early in the ontogeny. For instance, mutations in master regulatory genes start to exert their effects during early patterning of craniofacial development (e.g.
Tbx22 in mammals and zebrafish [47
], Endothelin1 in several vertebrates [48
]), while precocious treatment of Japanese flounder larvae (6-9 days post-hatching) with retinoid acid receptor agonists produces lower jaw deformities [49
]. In the present study, transcriptional changes are observed at a much later stage, which might represent either the downstream effects of earlier changes in the expression of hierarchically higher genes and/or transcriptional perturbations that start at a less differentiated stage and are maintained until later stages. On the other hand, no statistically significant differences in the expression profile of 38-days old sea bass juveniles were observed comparing deformed against normal fish. For this stage, however, whole heads were analyzed, which might have caused a reduced sensitivity in detecting differential expression. Using image analysis on juvenile sea bass pictures it was estimated that the lower jaw represents less than 1/10 of the entire head, while the anterior region of the mandible, the dentary, is approximately 1/30 of the sampled tissues. This means that if for instance a certain gene has an average expression level of 10 in the whole head whereas it is up-regulated 10-fold (100) exclusively in the mandibular region of deformed animals, the observed overall fold change would be 1.9, and in case differential expression is limited to the dentary a 1.3 fold increase would be present, likely below the threshold of reliable detection using microarray analysis. Microdissection of the mandible or part of it will therefore be required to identify differentially expressed genes in deformed animals at early stages, as reported here for 58-days old juveniles. The relatively large set of down-regulated transcripts in the lower jaw of prognathous fish also suggests that mandibular prognathism observed in the present study might be distinct from pugheadness, which is reported to be the result of under-development of the upper jaw, with only apparent protrusion of the normally developed lower jaw. As in most teleost fish, mandible development in D. labrax
is characterized by the presence of a cartilagenous component, the Meckel's cartilage, which is established early and subsequently regresses and ossifies, and a dermal bone component, the dentary, which appears later (25 days post hatching) and undergoes direct (membranous) ossification [50
]. Both components receive an important contribution of cells originating from the neural crests [51
]. The ossification process of the Meckel's cartilage is relatively well known, with a central role of the Hedgehog pathway, and it appears to be conserved in fish and mammals [52
], whereas the molecular mechanisms controlling the dentary development are less characterized and seem to be at least partially distinct from those observed in the Meckel's cartilage [53
]. How can this relate to the genes that were found to be differentially expressed in the present study? A significant enrichment in regulatory genes, especially in the development of anatomical structures was obtained in GO functional annotation, and gene-by-gene evaluation confirmed that some of them might be linked to bone formation. For instance, SOX4 has been reported to act downstream to Parathyroid Hormome (PTH) and PTH related protein (PTHrP) in osteoblast-like cells, being highly expressed in hypertrophic condrocytes during the mineralizing phase of endochondral ossification [54
]. Pleiotrophin or HB-GAM is highly expressed in bone, where it seems to play a role in bone development and remodelling [55
]. Retinoic acid receptor X gamma (RXRγ) forms heterodimers with vitamin D receptors, Peroxisome prolifireator-activated receptors (PPARs), and Retinoic acid receptors (RARs). All these nuclear receptors are somehow involved in bone formation and more in general in controlling skeletal growth. For instance, vitamin A and its metabolites, acting through RARs, have been reported to cause vertebral and craniofacial deformities in farmed fish [5
]. Remarkably, high levels of dietary vitamin A significantly increased the frequency of cranial malformations (underdeveloped lower jaw) in D. labrax
] as well as in other fish species [49
It seems, however, that bone development is not the most represented biological process when examining the functional roles of protein encoded by differentially expressed genes in deformed sea bass. A substantial fraction of them appears to affect nerve growth and function. This evidence might be explained by a temporal shift in the development of mandibular nerves, which in turn may be a consequence of the altered process of mandible bone formation. Down-regulation of several markers of neuronal development seems to suggest a delay in the differentiation of neuronal cells that constitute the nervous component of the lower jaw. However, there could be a closer relationship between mandibular bone formation and transcriptional changes for neuron-specific genes. Increasing evidence indicates that the nervous system participates in the regulation of bone physiology [57
]. Peripherally-released neurotransmitters exert their actions on osteoblasts and osteoclasts, which have been demonstrated to express specific receptors for these mediators. Of particular interest for the present study is the role of calcitonin gene-related peptide (CGRP), which is a neuropeptide with a well-established function in bone metabolism, promoting bone formation and repressing bone resorption [58
]. As other neuromediators, CGRP is also produced directly by osteoblasts as an autocrine factor [57
]. The expression of neuronal-specific genes in bone cells is not limited to neurotransmitters and their receptors, but it extends to the molecular network for regulated glutamate exocytosis (e.g.
SNARE, SNAP-25, syntaxin, synaptophysin, syntaghmin) generally described in pre-synaptic nerve terminals [59
]. The role of glutamate signalling in osteoblasts and osteoclasts is complex and in vivo
studies are still limited [60
]. Evidence from conditional KO mice lacking components of the glutamate pathways showed reduced mineralization and delayed ossification. Several genes involved in exocytosis and glutamate signalling are represented among down-regulated transcripts in the present study (e.g.
SNAP-25, Synaptophysin, Synapsin, Metabotropic glutamate receptor 8). A working hypothesis to hold together the above evidence might be that in jaw-protruding sea bass the bone formation process is delayed compared to normally developed animals, through down-regulation of different signalling pathways, which control bone formation/remodelling, either directly in osteoblasts and/or in neuronal terminals innervating the mandible. A delay in the ossification process might allow a prolonged growth of the mandibular bone components resulting in a protruding lower jaw.
Finally, it should be noted that part of the transcriptional differences observed in deformed animals might point, at least in part, to other regions of the dissected mandible. For instance, the formation of the tongue and/or the teeth might be indirectly affected by the deformity, and contribute to the differential gene profile. This in turn might provide a complementary/alternative hypothesis to explain the prevalence of neuron-related transcripts. In fact, tooth development has been shown to be tightly linked with nerve development.
Clearly, the above hypotheses await further confirmation from additional developmental stages and with the use of methods that allow better dissection of anatomical structures.