Gene expression profiling identifies different sub-types of retinoblastoma.

BACKGROUND: Mutation of the RB1 gene is necessary but not sufficient for the development of retinoblastoma. The nature of events occurring subsequent to RB1 mutation is unclear, as is the retinal cell-of-origin of this tumour.
METHODS: Gene expression profiling of 21 retinoblastomas was carried out to identify genetic events that contribute to tumorigenesis and to obtain information about tumour histogenesis.
RESULTS: Expression analysis showed a clear separation of retinoblastomas into two groups. Group 1 retinoblastomas express genes associated with a range of different retinal cell types, suggesting derivation from a retinal progenitor cell type. Recurrent chromosomal alterations typical of retinoblastoma, for example, chromosome 1q and 6p gain and 16q loss were also a feature of this group, and clinically they were characterised by an invasive pattern of tumour growth. In contrast, group 2 retinoblastomas were found to retain many characteristics of cone photoreceptor cells and appear to exploit the high metabolic capacity of this cell type in order to promote tumour proliferation.
CONCLUSION: Retinoblastoma is a heterogeneous tumour with variable biology and clinical characteristics.

Disruption of the RB pathway as a consequence of RB1 gene mutation or mutation of other pathway components (e.g., D-type cyclins, CDK4 or p16INK4A) is common to most if not all human cancers. Although these mutations are assumed to impact primarily on cell cycle regulation, RB1 has many additional roles, including regulation of chromosome stability, senescence and cellular differentiation (Indovina ). Regulation of permanent cell cycle withdrawal and terminal differentiation are especially relevant in a developmental context and may be important in the prototype RB1-associated tumour, retinoblastoma. An understanding of the molecular pathogenesis of retinoblastoma may also shed light on additional roles of pRB in other tumour types which also show high frequencies of RB1 gene mutation/deletion (rather than cyclin/CDK/CDKI mutations), for example, osteosarcoma, small-cell lung cancer, bladder cancer and hepatocellular carcinoma.

Knowledge of the cell-of-origin of retinoblastoma is important to understand the role of RB1 mutation in this tumour. The mammalian retina is composed of six neuronal cell types, rod, cone, horizontal, amacrine, bipolar and ganglion cells, and one glial cell type, Müller glia (Figure 1), all of which are derived from a common retinal progenitor cell (RPC) and arise in an evolutionarily conserved birth order during development. Retinal progenitor cells at a specific stage of development show competence to produce post-mitotic precursor cells, with the potential for terminal differentiation into a limited subset of retinal cell types, for example, early progenitor cells give rise to ganglion and cone precusor cells, while late progenitor cells give rise to Müller glial and bipolar precursor cells (Livesey and Cepko, 2001; Dyer and Bremner, 2005).

Figure 1

The human retina. (A) Organisation of the retina. (B) A H&E-stained section of the retina shows the outer and inner segments (OS, IS) of rod and cone photoreceptor cells. Photoreceptor nuclei form the outer nuclear layer (ONL). Nuclei of bipolar, amacrine, horizontal and Müller glial cells form the inner nuclear layer (INL), and the nuclei of ganglion cells form the ganglion cell layer (GCL). The outer plexiform layer (OPL) contains the processes and synaptic terminals of photoreceptors, horizontal and bipolar cells. The inner plexiform layer (IPL) contains the processes and synaptic terminals of bipolar, amacrine and ganglion cells. The processes of Müller glial cells fill all the retinal space not occupied by neurons and blood vessels. (Reproduced from Sung and Chuang (2010)). Abbreviations: A, amacrine cell; B, bipolar cell; C, cone; G, ganglion cell; H, horizontal cell; M, Müller glial cell; R, rod; RPE, retinal pigment epithelium.

Several studies of both human retinoblastoma and mouse retinoblastoma models have sought to define the retinal cell-of-origin of retinoblastoma and to understand why these cells are so susceptible to oncogenic transformation following RB1 mutation. These studies have reached a variety of conclusions, suggesting that retinoblastomas may be heterogeneous in their origin. A detailed immunohistochemical and genetic examination of human retinoblastoma led Xu to suggest, for example, that human retinoblastoma has properties of a cone precursor cell. Investigation of mouse models of retinoblastoma have generally implicated other cell types, however, including horizontal, amacrine and Müller glial precursors (Chen ; Ajioka ; Pajovic ). A further study concluded that multiple cell type-specific developmental programs, including those of amacrine, horizontal, photoreceptor and progenitor cells, are co-expressed in individual human retinoblastoma cells (McEvoy ).

A general finding from studies of both human and murine retinoblastoma is that secondary mutations in addition to loss of RB1 are required for tumour development. It has been reported, for example, that loss of RB1 in humans induced non-proliferative retinoma and that increasing genome instability was correlated with progression to retinoblastoma (Dimaras ). The observation of recurrent chromosomal alterations including 1q and 6p gain and 16q loss in up to 50% of all retinoblastomas suggests that genes on these chromosomes may contribute to tumour development and progression. A number of genes of interest have been highlighted, for example, MDM4 and KIF14 on 1q, E2F3 on 6p, CDH11 and RBL2/p130 on 16q (Corson and Gallie, 2007). However, the significance of many additional genes that show altered patterns of expression in association with these chromosomal alterations (e.g., CENPF, DEK, KIFC1) remains to be determined (Corson and Gallie, 2007).

In this study, we have used gene expression profiling to gain further insight into the molecular pathways which drive retinoblastoma tumorigenesis and their relationship to retinal developmental processes.

Materials and methods

Patient samples

Frozen tissue from 21 retinoblastomas, enucleated without prior treatment, was used for microarray analyses and qRT-PCR. Tissue was obtained from dissected globes and was divided for RNA and DNA extraction. Histopathology was assessed (blinded to genetic data) on tumour tissue, which was formalin-fixed and paraffin-embedded at the time of diagnosis. Five to 10 tissue sections were assessed according to the guidelines of the International Retinoblastoma Staging Working Group (Sastre ). Patients were all treated at the Birmingham's Children Hospital, UK, and consent was obtained for tissue banking for ethically approved research. The study was approved by the local Research Ethics Committee.

Array-CGH

DNA was prepared from retinoblastomas and peripheral blood lymphocytes from 10 normal healthy individuals using phenol/chloroform extraction (normal DNAs were then pooled for further processing). DNA was sonicated (100–1000 bp), purified (Qiagen Ltd, Crawley, UK) and then amplified in a two part procedure using primer A (5′-GTTTCCCAGTCACGGTCNNNNNNNNN-3′) for initial linear amplification, followed by primer B (5′-GTTTCCCAGTCACGGTC-3′) (Wang ) for 30 cycles of amplification with the inclusion of aminoallyl-dUTP to faciliate post amplification labelling. Fragmentation and labelling of 7.5 μg of amplified DNA was performed using the GeneChip WT dsDNA Terminal Labeling kit (Affymetrix UK Ltd, High Wycombe, UK). Hybridization to Affymetrix GeneChip Human Promoter 1.0R tiling arrays was carried out according to the standard Affymetrix procedures.

Gene expression arrays

Total RNA was isolated from retinoblastoma samples using Trizol extraction (Life Technologies Ltd, Paisley, UK), and additional purification was carried out using the Qiagen RNeasy system (Qiagen). RNA integrity was assessed using the Agilent Bioanalyzer 2100 (Agilent Technologies UK Ltd, Wokingham, UK). First- and second-strand cDNA synthesis, labelling and hybridisation to Affymetrix Human Gene 1.0ST arrays were performed according to the standard Affymetrix protocol. Affymetrix data were extracted, normalised and summarised using the robust multi-average method implemented in the Affymetrix Expression Console. CEL files from the same Affymetrix array platform (Human Gene 1.0ST) were downloaded from the GEO database (http://www.ncbi.nlm.nih.gov/geo/) for normal adult retina (GSM607947, GSM607948) and 96-day human fetal retina (GSM460264), and were normalised and summarised in parallel with retinoblastoma samples.

Array data analysis

CEL files from Affymetrix GeneChip Human Promoter 1.0R Tiling Arrays were analysed for copy number alterations using Partek Genomic Suite (PGS). Data were first log2 transformed and quantile normalised (with adjustment for GC content and probe sequence). Using the normal sample as baseline, copy number estimates were then obtained for each of 4.2 million probes on the array and segmentation applied with parameters set at: minimum number of genomic markers: 10, P-value threshold: 0.01 and signal to noise: 1.0.

Principal component and hierarchical clustering (HC) analysis of expression data from Affymetrix Human Gene 1.0ST arrays was carried out using PGS. Default parameters were used for principal component analysis (PCA; dispersion matrix: correlation; eigenvector scaling: normalised) and elipsoids were plotted at 2 standard deviations from the centroid of each group. Hierarchical clustering metrics used were Pearson correlation (distance measurement) and average linkage clustering. For other analyses, summarised Affymetrix data were imported into MultiExperiment Viewer (Saeed ). Differentially expressed genes were identified using SAM (Significance Analysis of Microarrays) with a q-value or false discovery cutoff at 5%. Functional annotation of lists of differentially expressed genes was carried out using the Database for Annotation, Visualisation and Integrated Discovery (DAVID) (Huang ). Annotation categories searched included GOTERM_BP (Biological Process), GOTERM_MF (Molecular Function) and pathways from the Biocarta, Kegg and Reactome databases. Hierarchical clustering of the reduced data set (retinal-associated genes) was performed using the GenePattern Software Suite (Reich ). Geneset enrichment analysis (GSEA) was carried out using canonical pathways included in the Molecular Signatures Database (v3.0) of the GSEA software package (Subramanian ). Data for GSEA was preprocessed so that only a single probeset for each gene (that with the highest expression) was included.

Chromosome 1q, 6p and 16q differentially expressed genes (fold change ⩾1.5, q⩽0.05) were defined as potentially cancer related on the basis of encoding proteins of known function and citations in the literature, implicating involvement in cancer or in cellular regulatory processes.

Real-time RT-PCR

One microgram of total RNA was used for cDNA synthesis using the QuantiTect Reverse Transcription kit (Qiagen) according to the manufacturer's protocol. Quantitative RT-PCR was performed in 1 × SYBR Green PCR mixture (Applied Biosystems, Warrington, UK), 1 μl diluted cDNA (1 : 5) and 0.5 μℳ gene-specific primers (Supplementary Table S7) designed using Primer 3 (http://frodo.wi.mit.edu) software. The reactions were performed in triplicate and β-actin was used as endogenous normalisation control.

Results

Retinoblastoma is a heterogeneous tumour

Initial exploratory analysis of the gene expression profiles of 21 retinoblastomas (Table 1) was carried out to assess heterogeneity. Principal component analysis showed a clear separation between normal retina and retinoblastoma, and in addition indicated a further sub-division of retinoblastomas into two main groups with 13 and 6 members, respectively (Figure 2A). A further two samples (group 3) were subsequently shown to be more similar to normal retina (see below). This grouping was also shown in a HC analysis (Figure 2B).

Table 1

Characteristics of retinoblastomas

										Chromosome copy number alterationsf
RB ID	Histological differentiation	Types of rosettes	Extent of rosettesa	Apoptosisb	Choroid invasionc	Choroid invasion depth	Optic nerve invasiond	Scleral invasion	Genetic groupe	1q gain	6p gain	16q loss
RB1	Poor	FW, N	1	2	2	D	3	0	1	+	−	+
RB3	Poor	N	1	2	0	S	1	0	1	+	+	−
RB4	Poor	N	2	2	1	D	3	0	1	−	+	+
RB6	Well	FW	3	2	0	S	0	0	1	−	(+)g	−
RB7	Poor	−	0	2	2	D	3	0	1	+	−	+
RB9	Poor	−	0	2	2	D	4	0	1	−	−	−
RB10	Well	FW, N	4	2	0	S	0	0	1	−	+	−
RB11	Poor	FW, P	1	2	0	S	2	0	1	−	+	−
RB16	Intermed	FW, N, P	2	2	2	D	3	0	1	ND	ND	ND
RB17	Poor	N	1	2	1	D	3	0	1	−	−	−
RB18	Poor	P	0	2	2	D	1	0	1	+	+	+
RB19	Poor	FW, N, P	1	2	2	D	1	0	1	−	+	−
RB20	Well	FW	3	2	0	S	0	0	1	−	+	−
RB2	Well	FW	4	1	1	D	3	0	2	−	−	−
RB5	Well	FW, N	3	1	0	S	0	0	2	−	−	−
RB8	Well	FW	3	1	2	D	1	1	2	−	−	−
RB14	Well	FW, N	4	1	0	S	0	0	2	ND	ND	ND
RB15	Intermed	FW	2	2	0	S	0	0	2	−	−	−
RB21	poor	FW, N, P	1	2	0	S	0	0	2	−	−	−
RB12	Well	FW, N	3	1	1	D	2	0	3	−	−	−
RB13	Well	FW, N	4	1	0	S	0	0	3	−	−	−

Abbreviations: FW=Flexner-Wintersteiner (characteristic of retinoblastoma); N=Homer-Wright rosettes (associated with tumours of neural origin); ND=not determined; P=pallisading; S, D=superficial, deep choroidal invasion.

Extent of rosettes: 0, none; 1, <10% 2, 10–25% 3, 25–50% 4, >50%.

Apoptosis: 1, low (patchy distribution); 2, high (observed over whole tumour area).

Choroid invasion: 0, none; 1, localised; 2, extensive.

Optic nerve invasion: 0, none; 1, pre-laminar; 2, intra-laminar; 3, retro-laminar; 4, to resection margin.

Genetic group–see text.

+/−: Presence/absence of chromosome copy number alteration; nd: not done.

Gain restricted to ∼3 Mb within 6p21.3.

Figure 2

Retina and retinoblastoma gene expression data. (A) Principal component analysis shows separation of retinoblastoma samples into three groups, group 1 (blue symbols), group 2 (red symbols) and group 3 (green symbols). Normal adult retina and fetal retina are from http://www.ncbi.nlm.nih.gov/geo/. (B) Hierarchical clustering (HC) of the same samples. Colour coding shows concordance with PCA.

Functional annotation of differentially expressed genes in retina and retinoblastma (fold change ⩾2.0, q⩽0.05; Supplementary Table S1) indicated that, overall, the most highly enriched gene ontology (GO) categories were those associated with photoreceptor differentiation/maintenance and phototransduction in normal adult retina, and with cell cycle and mitosis in retinoblastoma (Table 2). More detailed examination of the differences between the three groups, using pairwise comparisons (Supplementary Tables S2–S4), identified significant enrichment of GO categories relating to mitosis, spindle function, DNA replication and cytokinesis in groups 1 and 2 relative to group 3 (Table 3). Group 2 was additionally characterised by GO categories relating to photoreceptor development and differentiation, while in both groups 2 and 3 phototransduction was also highlighted (Table 4).

Table 2

Gene ontology of genes differentially expressed (a) in normal adult retina and (b) retinoblastoma

Cluster	GO term	Fold enrichment	P-value
(a)
Annotation cluster 1	Enrichment score: 14.46
GOTERM_BP_FAT	GO:0007601∼visual perception	3.32	0.0000
GOTERM_BP_FAT	GO:0050953∼sensory perception of light stimulus	3.32	0.0000
GOTERM_BP_FAT	GO:0007600∼sensory perception	1.07	0.3013
Annotation cluster 2	Enrichment score: 8.82
GOTERM_BP_FAT	GO:0008104∼protein localisation	1.65	0.0000
GOTERM_BP_FAT	GO:0045184∼establishment of protein localisation	1.65	0.0000
GOTERM_BP_FAT	GO:0015031∼protein transport	1.64	0.0000
Annotation cluster 3	Enrichment score: 7.20
GOTERM_BP_FAT	GO:0007602∼phototransduction	5.64	0.0000
GOTERM_BP_FAT	GO:0009583∼detection of light stimulus	5.01	0.0000
GOTERM_BP_FAT	GO:0009582∼detection of abiotic stimulus	3.29	0.0000
GOTERM_BP_FAT	GO:0009581∼detection of external stimulus	2.89	0.0000
GOTERM_BP_FAT	GO:0051606∼detection of stimulus	2.37	0.0000
Annotation cluster 4	Enrichment score: 6.15
GOTERM_BP_FAT	GO:0046530∼photoreceptor cell differentiation	5.58	0.0000
GOTERM_BP_FAT	GO:0042461∼photoreceptor cell development	5.43	0.0000
GOTERM_BP_FAT	GO:0001754∼eye photoreceptor cell differentiation	5.32	0.0000
GOTERM_BP_FAT	GO:0042462∼eye photoreceptor cell development	5.12	0.0000
Annotation cluster 5	Enrichment score: 5.25
GOTERM_MF_FAT	GO:0043167∼ion binding
GOTERM_MF_FAT	GO:0046872∼metal ion binding	1.17	0.0000
GOTERM_MF_FAT	GO:0043169∼cation binding	1.16	0.0000
GOTERM_MF_FAT	GO:0043169∼cation binding	1.16	0.0000
(b)
Annotation cluster 1	Enrichment score: 185.65
GOTERM_MF_FAT	GO:0004984∼olfactory receptor activity	9.23	0.0000
GOTERM_BP_FAT	GO:0007606∼sensory perception of chemical stimulus	7.95	0.0000
GOTERM_BP_FAT	GO:0007608∼sensory perception of smell	8.35	0.0000
KEGG_PATHWAY	hsa04740:Olfactory transduction	6.36	0.0000
GOTERM_BP_FAT	GO:0007600∼sensory perception	4.96	0.0000
GOTERM_BP_FAT	GO:0007186∼G-protein-coupled receptor protein signalling pathway	4.09	0.0000
GOTERM_BP_FAT	GO:0050890∼cognition	4.46	0.0000
GOTERM_BP_FAT	GO:0050877∼neurological system process	3.50	0.0000
Annotation cluster 2	Enrichment score: 20.53
GOTERM_BP_FAT	GO:0000279∼M phase	3.44	0.0000
GOTERM_BP_FAT	GO:0022403∼cell cycle phase	2.94	0.0000
GOTERM_BP_FAT	GO:0022402∼cell cycle process	2.30	0.0000
GOTERM_BP_FAT	GO:0007049∼cell cycle	2.07	0.0000
Annotation cluster 3	Enrichment score: 19.54
GOTERM_BP_FAT	GO:0000280∼nuclear division	3.68	0.0000
GOTERM_BP_FAT	GO:0007067∼mitosis	3.68	0.0000
GOTERM_BP_FAT	GO:0000087∼M phase of mitotic cell cycle	3.62	0.0000
GOTERM_BP_FAT	GO:0048285∼organelle fission	3.54	0.0000
GOTERM_BP_FAT	GO:0000278∼mitotic cell cycle	2.64	0.0000
Annotation cluster 4	Enrichment score: 7.86
GOTERM_BP_FAT	GO:0051321∼meiotic cell cycle	3.46	0.0000
GOTERM_BP_FAT	GO:0051327∼M phase of meiotic cell cycle	3.41	0.0000
GOTERM_BP_FAT	GO:0007126∼meiosis	3.41	0.0000
Annotation cluster 5	Enrichment score: 5.79
GOTERM_BP_FAT	GO:0065004∼protein-DNA complex assembly	3.67	0.0000
GOTERM_BP_FAT	GO:0031497∼chromatin assembly	3.56	0.0000
GOTERM_BP_FAT	GO:0006334∼nucleosome assembly	3.41	0.0000
GOTERM_BP_FAT	GO:0034728∼nucleosome organisation	3.20	0.0000
GOTERM_BP_FAT	GO:0006333∼chromatin assembly or disassembly	2.63	0.0000
GOTERM_BP_FAT	GO:0006325∼chromatin organisation	1.07	0.4285

Abbreviation: GO=gene ontology.

Genes with at least two-fold (q⩽0.05) differential expression in normal adult retina vs retinoblastoma were analysed using DAVID (http://david.abcc.ncifcrf.gov). The top five most highly enriched clusters of GO terms are presented. The group enrichment score used to rank biological significance is the geometric mean (in −log scale) of member's P-values in a corresponding annotation cluster. GOTERM_BP: Biological Process. GOTERM_MF: Molecular function. G0_FAT: includes a subset of the Gene Ontology term set created in order to filter the broadest terms so that they do not overshadow the more specific terms.

Table 3

Gene ontology of genes differentially expressed in (a) group 1 and (b) group 2, both relative to group 3

Cluster	GO Term	Fold enrichment	P-value
(a)
Annotation cluster 1	Enrichment score: 34.87
GOTERM_BP_FAT	GO:0000279∼M phase	8.99	0.0000
GOTERM_BP_FAT	GO:0022403∼cell cycle phase	7.49	0.0000
GOTERM_BP_FAT	GO:0007067∼mitosis	10.46	0.0000
GOTERM_BP_FAT	GO:0000280∼nuclear division	10.46	0.0000
GOTERM_BP_FAT	GO:0000087∼M phase of mitotic cell cycle	10.28	0.0000
GOTERM_BP_FAT	GO:0048285∼organelle fission	10.05	0.0000
GOTERM_BP_FAT	GO:0022402∼cell cycle process	5.82	0.0000
GOTERM_BP_FAT	GO:0000278∼mitotic cell cycle	7.11	0.0000
Annotation cluster 2	Enrichment score: 29.41
GOTERM_MF_FAT	GO:0004984∼olfactory receptor activity	9.12	0.0000
GOTERM_BP_FAT	GO:0007608∼sensory perception of smell	7.74	0.0000
GOTERM_BP_FAT	GO:0007606∼sensory perception of chemical stimulus	7.17	0.0000
KEGG_PATHWAY	hsa04740:Olfactory transduction	5.66	0.0000
GOTERM_BP_FAT	GO:0007600∼sensory perception	4.29	0.0000
GOTERM_BP_FAT	GO:0050890∼cognition	3.82	0.0000
GOTERM_BP_FAT	GO:0007186∼G-protein-coupled receptor protein signalling pathway	3.39	0.0000
GOTERM_BP_FAT	GO:0050877∼neurological system process	2.99	0.0000
GOTERM_BP_FAT	GO:0007166∼cell surface receptor linked signal transduction	2.15	0.0000
Annotation cluster 3	Enrichment score: 10.49
GOTERM_BP_FAT	GO:0065004∼protein-DNA complex assembly	11.36	0.0000
GOTERM_BP_FAT	GO:0031497∼chromatin assembly	10.80	0.0000
GOTERM_BP_FAT	GO:0006334∼nucleosome assembly	10.62	0.0000
GOTERM_BP_FAT	GO:0034728∼nucleosome organisation	9.60	0.0000
GOTERM_BP_FAT	GO:0006333∼chromatin assembly or disassembly	7.40	0.0000
GOTERM_BP_FAT	GO:0034621∼cellular macromolecular complex subunit organisation	3.55	0.0000
GOTERM_BP_FAT	GO:0034622∼cellular macromolecular complex assembly	3.69	0.0000
GOTERM_BP_FAT	GO:0006325∼chromatin organisation	2.73	0.0001
Annotation cluster 4	Enrichment score: 8.65
GOTERM_BP_FAT	GO:0051327∼M phase of meiotic cell cycle	7.67	0.0000
GOTERM_BP_FAT	GO:0007126∼meiosis	7.67	0.0000
GOTERM_BP_FAT	GO:0051321∼meiotic cell cycle	7.52	0.0000
Annotation CLUSTER 5	Enrichment score: 8.27
GOTERM_BP_FAT	GO:0006974∼response to DNA damage stimulus	3.90	0.0000
GOTERM_BP_FAT	GO:0006281∼DNA repair	4.47	0.0000
GOTERM_BP_FAT	GO:0033554∼cellular response to stress	2.66	0.0000
(b)
Annotation cluster 1	Enrichment score: 66.78
GOTERM_BP_FAT	GO:0022403∼cell cycle phase	10.48	0.0000
GOTERM_BP_FAT	GO:0000279∼M phase	11.91	0.0000
GOTERM_BP_FAT	GO:0022402∼cell cycle process	8.34	0.0000
Annotation cluster 2	Enrichment score: 53.36
GOTERM_BP_FAT	GO:0000278∼mitotic cell cycle	9.83	0.0000
GOTERM_BP_FAT	GO:0007067∼mitosis	13.57	0.0000
GOTERM_BP_FAT	GO:0000280∼nuclear division	13.57	0.0000
GOTERM_BP_FAT	GO:0000087∼M phase of mitotic cell cycle	13.33	0.0000
GOTERM_BP_FAT	GO:0048285∼organelle fission	13.04	0.0000
Annotation cluster 3	Enrichment score: 21.83
GOTERM_BP_FAT	GO:0006281∼DNA repair	7.23	0.0000
GOTERM_BP_FAT	GO:0006974∼response to DNA damage stimulus	6.13	0.0000
GOTERM_BP_FAT	GO:0033554∼cellular response to stress	4.29	0.0000
Annotation cluster 4	Enrichment score: 14.86
GOTERM_BP_FAT	GO:0007126∼meiosis	10.47	0.0000
GOTERM_BP_FAT	GO:0051327∼M phase of meiotic cell cycle	10.47	0.0000
GOTERM_BP_FAT	GO:0051321∼meiotic cell cycle	10.26	0.0000
Annotation cluster 5	Enrichment score: 13.93
GOTERM_BP_FAT	GO:0065004∼protein-DNA complex assembly	12.82	0.0000
GOTERM_BP_FAT	GO:0031497∼chromatin assembly	12.33	0.0000
GOTERM_BP_FAT	GO:0034728∼nucleosome organisation	11.04	0.0000
GOTERM_BP_FAT	GO:0006334∼nucleosome assembly	11.66	0.0000
GOTERM_BP_FAT	GO:0006333∼chromatin assembly or disassembly	8.45	0.0000
GOTERM_BP_FAT	GO:0034622∼cellular macromolecular complex assembly	4.40	0.0000
GOTERM_BP_FAT	GO:0006325∼chromatin organisation	3.21	0.0000

Abbreviation: GO=gene ontology.

Genes with at least two-fold (q⩽0.05) differential expression in group 1 or group 2 vs group 3 were analysed using DAVID (http://david.abcc.ncifcrf.gov). The top five most highly enriched clusters of GO terms are presented.

Table 4

Gene ontology of genes differentially expressed in (a) group 2 and (b) group 3 relative to group 1

Cluster	GO term	Fold enrichment	P-value
(a)
Annotation cluster 1	Enrichment score: 9.60
GOTERM_BP_FAT	GO:0050953∼sensory perception of light stimulus	8.45	0.0000
GOTERM_BP_FAT	GO:0007601∼visual perception	8.45	0.0000
GOTERM_BP_FAT	GO:0050890∼cognition	2.28	0.0000
GOTERM_BP_FAT	GO:0007600∼sensory perception	2.33	0.0000
GOTERM_BP_FAT	GO:0050877∼neurological system process	1.98	0.0001
Annotation cluster 2	Enrichment score: 5.52
GOTERM_BP_FAT	GO:0007602∼phototransduction	19.07	0.0000
GOTERM_BP_FAT	GO:0009583∼detection of light stimulus	16.13	0.0000
GOTERM_BP_FAT	GO:0009582∼detection of abiotic stimulus	9.68	0.0000
GOTERM_BP_FAT	GO:0009581∼detection of external stimulus	8.50	0.0000
GOTERM_BP_FAT	GO:0009416∼response to light stimulus	5.47	0.0000
GOTERM_BP_FAT	GO:0051606∼detection of stimulus	5.33	0.0001
GOTERM_BP_FAT	GO:0009314∼response to radiation	3.78	0.0003
GOTERM_BP_FAT	GO:0009628∼response to abiotic stimulus	2.74	0.0008
Annotation cluster 3	Enrichment score: 4.19
GOTERM_BP_FAT	GO:0042462∼eye photoreceptor cell development	18.88	0.0000
GOTERM_BP_FAT	GO:0001754∼eye photoreceptor cell differentiation	17.98	0.0000
GOTERM_BP_FAT	GO:0042461∼photoreceptor cell development	15.73	0.0000
GOTERM_BP_FAT	GO:0046530∼photoreceptor cell differentiation	15.10	0.0000
GOTERM_BP_FAT	GO:0048592∼eye morphogenesis	5.47	0.0047
Annotation cluster 4	Enrichment score: 3.06
GOTERM_BP_FAT	GO:0048858∼cell projection morphogenesis	3.60	0.0001
GOTERM_BP_FAT	GO:0032990∼cell part morphogenesis	3.44	0.0002
GOTERM_BP_FAT	GO:0032989∼cellular component morphogenesis	2.38	0.0043
GOTERM_BP_FAT	GO:0000902∼cell morphogenesis	2.47	0.0044
Annotation cluster 5	Enrichment score: 1.96
GOTERM_BP_FAT	GO:0010324∼membrane invagination	2.86	0.0084
GOTERM_BP_FAT	GO:0006897∼endocytosis	2.86	0.0084
GOTERM_BP_FAT	GO:0016044∼membrane organisation	2.15	0.0183
(b)
Annotation cluster 1	Enrichment score: 20.08
GOTERM_BP_FAT	GO:0050953∼sensory perception of light stimulus	5.17	0.0000
GOTERM_BP_FAT	GO:0007601∼visual perception	5.17	0.0000
GOTERM_BP_FAT	GO:0007600∼sensory perception	1.59	0.0001
Annotation cluster 2	Enrichment score: 10.44
GOTERM_BP_FAT	GO:0007602∼phototransduction	10.40	0.0000
GOTERM_BP_FAT	GO:0009583∼detection of light stimulus	8.80	0.0000
GOTERM_BP_FAT	GO:0009582∼detection of abiotic stimulus	5.81	0.0000
GOTERM_BP_FAT	GO:0009581∼detection of external stimulus	5.10	0.0000
GOTERM_BP_FAT	GO:0051606∼detection of stimulus	3.78	0.0000
GOTERM_BP_FAT	GO:0009416∼response to light stimulus	3.36	0.0000
Annotation cluster 3	Enrichment score: 6.40
GOTERM_BP_FAT	GO:0031644∼regulation of neurological system process	3.14	0.0000
GOTERM_BP_FAT	GO:0051969∼regulation of transmission of nerve impulse	3.15	0.0000
GOTERM_BP_FAT	GO:0050804∼regulation of synaptic transmission	3.16	0.0000
Annotation cluster 4	Enrichment score: 6.03
GOTERM_BP_FAT	GO:0048858∼cell projection morphogenesis	2.73	0.0000
GOTERM_BP_FAT	GO:0032990∼cell part morphogenesis	2.62	0.0000
GOTERM_BP_FAT	GO:0030030∼cell projection organisation	2.24	0.0000
GOTERM_BP_FAT	GO:0048667∼cell morphogenesis involved in neuron differentiation	2.63	0.0000
GOTERM_BP_FAT	GO:0031175∼neuron projection development	2.41	0.0000
GOTERM_BP_FAT	GO:0048812∼neuron projection morphogenesis	2.50	0.0000
GOTERM_BP_FAT	GO:0007409∼axonogenesis	2.40	0.0001
Annotation cluster 5	Enrichment score: 4.88
GOTERM_BP_FAT	GO:0008104∼protein localisation	1.67	0.0000
GOTERM_BP_FAT	GO:0045184∼establishment of protein localisation	1.65	0.0000
GOTERM_BP_FAT	GO:0015031∼protein transport	1.62	0.0000

Abbreviation: GO=gene ontology.

Genes with at least two-fold (q⩽0.05) differential expression in group 2 or group 3 vs group 1 were analysed using DAVID (http://david.abcc.ncifcrf.gov). The top five most highly enriched clusters of GO terms are presented.

The conclusion drawn from these results is that groups 1 and 2 appear to have the greatest proliferative potential, but these groups differ in the extent of photoreceptor differentiation. The expression profile of group 2 retinoblastomas appears to be most consistent with derivation from a retinal cone photoreceptor lineage. Group 1 retinoblastomas may derive either from a different retinal lineage or from an early uncommitted cell type, for example, a RPC.

Retinoblastomas show different patterns of retinal gene expression

To investigate the histogenesis of retinoblastoma groups in more detail, data were extracted for all genes on the array associated with retinal development and function. Genes were identified using the Affymetrix gene description; also included were genes listed by Hennig ), which function in the regulation of photoreceptor gene expression and genes listed by Byerly and Blackshaw (2009) with roles in vertebrate retinal development. Unsupervised HC of tumours, based on 80 retina-associated genes, produced an identical grouping to that obtained with PCA and with HC of the complete data set (Figure 3), and provided information about the expression of genes associated with different retinal cell types in each of the tumour groups.

Figure 3

Hierarchical clustering of retina and retinoblastoma based on 80 genes associated with retinal development and function. Abbreviations: aRET, adult retina; fRET, fetal retina. *Cone-enriched genes; **Rod-enriched genes.

The expression of cone photoreceptor-enriched transcription factors, RXRG and THRB, was of particular interest in view of the reported cone precursor origin of retinoblastoma (Xu ). The expression of RXRG and THRB was increased 1.4- and 1.5-fold, respectively, in group 2 retinoblastomas relative to group 1, and was three- to five-fold increased relative to normal adult retina. Furthermore, the high expression in group 2 of downstream genes encoding cone opsins (OPN1MW, OPN1LW), cone arrestin (ARR3) and cone phosphodiesterases (PDE6C, PDE6H) is consistent with the development of group 2 retinoblastomas from a cone lineage cell. The fact that there was little or no expression of either rod rhodopsin (RHO) or rod arrestin (SAG) in this group confirms that the observed gene expression profile is not a consequence of the inclusion of normal retinal tissue. Expression of RHO and SAG were 170- and 50-fold higher, respectively, in normal retina, where rods outnumber cones by ∼20 : 1. In contrast, the two group 3 samples (RB12 and RB13) showed high-level expression of both cone and rod genes (e.g., OPN1MW/LW, RHO, NR2E3), as well as several genes characteristic of other retinal neuronal and glial cell types (HES1, PAX6, PROX1, VSX2, RLBP1/CRALBP) (Livesey and Cepko, 2001; Cid ; Joly ). Their similarity to the normal retinal samples (Figure 3) is suggestive of inclusion of some normal retinal tissue, and these two samples were not analysed further.

Although the expression of RXRG and THRB was moderately reduced rather than absent in group 1, a much greater decrease (six- to eight-fold) in the expression of OPN1MW/LW genes relative to group 2 suggests that group 1 retinoblastomas do not complete the cone differentiation programme and may be arrested at an earlier stage of development. It was noted that some group 1 retinoblastomas also showed variable upregulation (⩽1.5-fold) of a number of genes encoding transcription factors implicated more generally in eye and retinal development, particularly in the development of retinal horizontal and amacrine cells, for example, FOXN4, PTFIA and MITF (Fujitani ). These results suggest that group 1 retinoblastomas have the potential for at least the early stages of differentiation of more than one retinal cell type.

To further validate these results, quantitative RT-PCR was carried out for selected genes characteristic of each of the three groups. In the case of all genes tested (FOXN4, ARR3, RXRG, RHO, NR2E3 and PAX6), expression patterns were in agreement with the array results (Supplementary Figure S1).

Histopathology and cytogenetics of retinoblastoma sub-types

We next examined whether the genetic classification of retinoblastomas into different groups, with differing levels of cone photoreceptor differentiation, was related to histopathological or cytogenetic characteristics. A detailed histopathological review of all retinoblastomas (carried out by a reviewer without knowledge of the genetic data) showed relatively good agreement between genetic and histological definitions of photoreceptor differentiation (Table 1). Four out of six (66.7%) group 2 retinoblastomas were classified as well differentiated on the basis of histological features (the extent of photoreceptor-like Flexner-Wintersteiner rosettes), compared with only 3 out of 13 (23.1%) group 1 retinoblastomas. However, these figures may under-estimate the extent of correlation between genetic and histopathological data, as tissue used for genetic analysis was representative of only a small part of the tumour. It is also of interest that group 1 retinoblastomas showed an increased frequency of both post-laminar optic nerve invasion (PLONI, 6 out of 13, 46.1%) and deep choroid invasion (DCI, 8 out of 13, 61.5%) compared with group 2 (1 out of 6, 16.7% and 2 out of 8, 25.0% for PLONI and DCI respectively), although these differences did not reach statistical significance (Table 1). Patients with these features are considered to be at higher risk of extra-ocular tumour spread (Shields ). Surprisingly, the frequency of apoptotic cells was also significantly higher in group 1 retinoblastomas (P=0.004).

The occurrence of cytogenetic alterations most frequently reported to be associated with adverse clinical features was also investigated, that is, gain of chromosomes 1q and 6p and loss of chromosome 16q. Array-CGH of 19 samples for which material was available showed that these alterations tended to occur together and were significantly more frequent in group 1 retinoblastomas (P=0.003 for any alteration; P=0.029 for 6p gain) (Table 1, Figure 4). Overall, these results lend further weight to the classification of retinoblastoma into two distinct groups with differing genetic and clinical characteristics.

Figure 4

Array-CGH of retinoblastomas. The chromosomes with the most frequent copy number alterations were chromosomes 1, 6 and 16. These alterations were restricted to group 1 retinoblastomas.

Gene expression changes associated with recurrent chromosome alterations

We next asked whether the altered pattern of differentiation in group 1 retinoblastomas might be related to the expression of genes on chromosomes 1q, 6p and 16q. SAM analyses identified a total of 555 genes mapping to these chromosomal regions, which showed significantly altered expression (fold change ⩾1.5; FDR q⩽0.05) in retinoblastomas with the relevant chromosomal alteration (Supplementary Table S5). A subset of 56 protein-coding genes had functions which potentially could contribute to retinoblastoma tumorigenesis (based on involvement in other tumour types or in cell regulatory processes) or which are relevant to retinal/photoreceptor development and function (Table 5).

Table 5

Differentially expressed genes on chromosomes 1q, 6p and 16q

Gene symbol	Gene name	Cyto band	Fold change	q (%)
CKS1B	CDC28 protein kinase regulatory subunit 1B	1q21.2	2.1	1.3
HDGF	Hepatoma-derived growth factor	1q23.1	1.7	0.0
VANGL2	Vang-like 2 (van gogh, Drosophila)	1q23.2	1.5	1.3
NUF2	NUF2, NDC80 kinetochore complex component	1q23.3	2.7	0.8
CENPL	Centromere protein L	1q25.1	1.6	0.9
CACYBP	Calcyclin-binding protein	1q25.1	1.6	0.7
ASPM	Asp (abnormal spindle) homolog,	1q31	2.6	2.0
KIF14	Kinesin family member 14	1q32.1	2.3	1.0
NR5A2	Nuclear receptor subfamily 5, group A, member 2	1q32.1	1.8	0.5
YOD1	YOD1 OTU deubiquinating enzyme 1	1q32.2	1.5	0.0
NEK2	NIMA (never in mitosis gene a)-related kinase 2	1q32.3	2.0	4.1
CENPF	Centromere protein F (mitosin)	1q41	2.4	2.7
PRIM2	Primase, DNA, polypeptide 2	6p11.2	1.7	0.4
MCM3	Minichromosome maintenance complex component 3	6p12	1.7	0.3
BMP5	Bone morphogenetic protein 5	6p12.1	1.7	0.0
TFAP2D	Transcription factor AP-2 delta	6p12.1	1.6	0.2
OPN5	Opsin 5	6p12.3	2.0	0.0
RUNX2	Runt-related transcription factor 2	6p21	1.7	0.0
UHRF1BP1	UHRF1-binding protein 1	6p21	1.5	0.0
CDC5L	CDC5 cell division cycle 5-like	6p21	1.5	0.8
STK38	Serine/threonine kinase 38	6p21	1.9	0.1
USP49	Ubiquitin-specific peptidase 49	6p21	1.8	0.1
POLH	Polymerase (DNA directed), eta	6p21.1	1.7	0.1
FOXP4	Forkhead box P4	6p21.1	1.5	0.1
PIM1	Pim-1 oncogene	6p21.2	1.7	0.0
TEAD3	TEA domain family member 3	6p21.2	1.6	0.0
DAXX	Death-domain associated protein	6p21.3	1.6	0.0
KIFC1	Kinesin family member C1	6p21.3	2.2	0.0
EHMT2	euchromatic histone-lysine N-methyltransferase 2	6p21.31	1.5	0.1
FANCE	Fanconi anaemia, complementation group E	6p21.31	1.8	0.0
E2F3	E2F transcription factor 3	6p22	1.5	0.2
JARID2	Jumonji, AT rich interactive domain 2	6p22.3	1.5	0.2
DEK	DEK oncogene	6p22.3	1.5	0.4
CDKAL1	CDK5 regulatory subunit associated protein 1-like 1	6p22.3	1.7	0.1
NUP153	Nucleoporin 153 kDa	6p22.3	1.6	0.1
SOX4	SRY-box 4	6p22.3	1.6	0.1
RNF182	Ring finger protein 182	6p23	1.5	0.2
PRPF4B	PRP4 pre-mRNA processing factor 4 homolog B	6p25.2	1.5	0.3
RNF8	RING finger protein 8	6p31.3	1.6	0.0
N4BP1	NEDD4-binding protein 1	16q12.1	2.0	1.2
CHD9	Chromodomain helicase DNA-binding protein 9	16q12.2	3.1	0.0
RBL2	retinoblastoma-like 2 (p130)	16q12.2	2.4	1.2
RPGRIP1L	RPGRIP1-like	16q12.2	2.3	2.4
RSPRY1	Ring finger and SPRY domain containing 1	16q13	2.2	0.8
BBS2	Bardet-Biedl syndrome 2	16q21	3.0	0.0
CNOT1	CCR4-NOT transcription complex, subunit 1	16q21	2.1	0.0
NAE1	NEDD8 activating enzyme E1	16q22	2.0	4.2
PRMT7	Protein arginine methyltransferase 7	16q22.1	1.6	1.2
SF3B3	Splicing factor 3b, subunit 3,	16q22.1	1.9	1.2
TERF2	Telomeric repeat-binding factor 2	16q22.1	1.9	4.1
WWP2	WW domain containing E3 ubiquitin protein ligase 2	16q22.1	1.7	2.4
TXNL4B	Thioredoxin-like 4B	16q22.2	1.9	4.1
CDK10	Cyclin-dependent kinase 10	16q24	1.5	0.8
USP10	Ubiquitin-specific peptidase 10	16q24.1	2.5	4.1
FBXO31	F-box protein 31	16q24.2	1.6	1.0
TCF25	Transcription factor 25 (basic helix-loop-helix)/NULP1	16q24.3	2.5	1.1

Chromosome 1q genes identified in retinoblastomas with 1q gain included several with functions in mitosis, particularly centrosome duplication, chromosome segregation and the mitotic spindle checkpoint (ASPM, CENPF, KIF14, NUF2, NEK2). Additional chromosome 1q genes identified have functions in promoting cell cycle progression (CKS1B) and transcriptional regulation during development (NR5A2).

Chromosome 6p genes that showed increased expression in association with 6p gain again emphasised the same functional categories highlighted by 1q genes, for example, mitosis (KIFC1, NUP153), cell cycle progression (E2F3, UHRF1BP1) and developmental transcriptional regulation (FOXP4, SOX4, TEAD3). TEAD3 is a downstream target of the Hippo pathway, and it is of interest that combined inactivation of RB and Hippo pathway kinases has been shown to result in dedifferentiation in the Drosophila retina (Nicolay ). Chromosome 6p genes also included several with roles in the epigenetic regulation of gene expression, particularly through modification of chromatin (DAXX, DEK, EHMT2, JARID2, RUNX2, PRPF4B).

Genes which were downregulated in association with 16q loss included RBL2 (p130), which has an essential role in maintaining permanent cell cycle arrest and senescence following events initiated by RB1 (Helmbold ). It is of interest that RBL2 expression was reduced not only in retinoblastomas with 16q loss (N=4), but also showed a 2.2-fold decrease overall in group 1 retinoblastomas relative to group 2, which may reflect epigenetic regulation (De Falco and Giordano, 2006). Other 16q genes identified included USP10, which is upregulated in response to DNA damage and which stabilises p53, also FBXO31 and NAE1, which function in ubiquitin pathways to regulate cell cycle progression, cell growth and survival, and CHD9 and PRMT7, which have roles in chromatin remodelling, with the latter also implicated in the DNA damage response (Karkhanis ). A further consequence of chromosome 16q loss was the decreased expression of BBS2 and RPGRIP1L, which are required for retinal photoreceptor function.

Geneset enrichment analysis analysis of regulatory pathways in group 1 and group 2 retinoblastomas

Geneset enrichment analysis was used to obtain further information about the regulatory pathways responsible for the differing characteristics of group 1 and 2 retinoblastomas. Genesets enriched in group 1 retinoblastomas included several related to G-protein-coupled receptor-mediated signalling (Supplementary Table S6). There was a strong emphasis on neurotransmitter signalling, particularly involving class A1 (RHO-like) GPCRs, which include the amine ligand-binding receptors (e.g., serotonin, adrenergic and cholinergic receptors) and peptide ligand-binding receptors (hypocretin/orexin, NPY, melanocortin and tachykinin receptors) (Supplementary Figure S2). These results are of interest as although glutamine is the major neurotransmitter within the retina, particularly for photoreceptor, bipolar and ganglion cells, additional neurotransmitters/receptors are a feature of other cell types. Serotonin and NPY, for example, are associated with subsets of amacrine cells (Bagnoli ; Ghai ), tachykinins with amacrine and ganglion cells, and orexins with horizontal, bipolar, amacrine and ganglion cells (Kolb ; Savaskan ). This genome-wide analysis of gene expression patterns is consistent with the idea that group 1 retinoblastomas show characteristics of multiple retinal cell types.

In the case of group 2 retinoblastomas, many significantly enriched genesets were relevant to photoreceptor function, particularly mitochondrial function/energy metabolism as well as transcription and RNA splicing (Supplementary Table S6, Supplementary Figure S3). Photoreceptors, especially cones, are among the most metabolically active cells in the body (Eckmiller 2004; Reidel ), and normal retinal/photoreceptor functions maintained in group 2 retinoblastomas may promote tumour growth.

Discussion

The aetiology of retinoblastoma and its relationship to retinal histogenesis have been intensively studied. However, as yet there is still relatively little consensus on the cell-of-origin or the molecular mechanisms involved. Our results suggest that retinoblastomas, like most other tumour types, are not a homogeneous group but may be classified into at least two different sub-types based on molecular and histopathological characteristics. Furthermore we suggest that these sub-types may have differing aetiology and origin.

Classification of retinoblastomas

Our results, based on gene expression profiling, are consistent with the existence of two sub-types of retinoblastoma, group 1 retinoblastomas that show features of multiple retinal cell types, and group 2 retinoblastomas that show a distinctive cone photoreceptor expression profile. This sub-division is further supported by the finding that recurrent chromosomal alterations typical of retinoblastoma (1q loss, 6p gain and 16q loss) were restricted to group 1 tumours, and in addition suggest that the genetic mechanisms driving tumorigenesis in each group may also be distinct. It was also noted that adverse histopathological features (PLONI and deep choroid invasion) were more frequent in group 1 retinoblastomas.

Our classification of retinoblastoma is in contrast to two previous reports, which emphasised the overall similarity between retinoblastomas, but nevertheless also provides support for both of these studies (Xu ; McEvoy ). Xu proposed a cone precursor origin for retinoblastoma (similar to our group 2 tumours), following immunohistochemical detection of cone-specific markers RXRG, THRB and cone arrestin in all 40 retinoblastomas examined. Consistent with this finding, both our group 1 and 2 retinoblastomas expressed genes encoding these proteins, albeit at higher levels in group 2. Similarly, the failure to detect several developmental markers of other retinal cell types (PROX1, CHX10/VSX2, PAX6) in retinoblastoma cells in the study of Xu et al (2009) is also consistent with our observation that these genes were significantly downregulated in both group 1 and 2 retinoblastomas relative to normal adult and normal fetal retina (Figure 3). We conclude, therefore, that the results of Xu et al (2009) do not exclude the possibility that some retinoblastomas also have characteristics of non-cone retinal cell types.

McEvoy concluded from gene expression analyses that retinoblastoma tumour cells express genes associated with multiple retinal cell types. This is similar to our findings for group 1 retinoblastomas. However, further examination of the data of McEvoy et al (2009; GEO:GSE29683) on 52 human retinoblastomas suggests that there is a significant range of expression values for individual retinal cell type-specific genes in this large complex data set (Supplementary Figure S4). These results are not inconsistent with the idea that some retinoblastomas have a predominantly cone-like expression profile.

The retinoblastoma cell-of-origin

The identification of different sub-types of retinoblastoma raises the question as to whether these represent divergent evolution from a common tumour progenitor cell or alternatively if there is a different cell-of-origin in each case. It is important to consider this question in the context of normal retinal development. Retinal progenitor cells progress through a series of competence states during development, giving rise to retinal cell type-specific precursors in a temporally restricted manner, a process which may be regulated at least in part by the sequential expression of retinal transcription factors (Ohsawa and Kageyama, 2008; Agathocleous and Harris, 2009). However, there is considerable overlap in the generation of different cell types and individual progenitor cells from the same developmental time point not only show extensive heterogeneity of gene expression but may also express genes encoding multiple retinal transcription factors within a single cell (Trimarchi ). Consistent with these observations, recent evidence suggests that RPCs show a surprising level of stochasticity in cell fate choice, and the observed variability in the number and type of cells resulting from clonal expansion of individual RPCs (a RPC may divide to produce two daughter RPCs (P/P), a RPC and a differentiating cell (P/D) or two differentiating cells (not necessarily identical) (D/D)) fits well with a model in which cells divide or differentiate with fixed probabilities related to the ratios of different cell types in the normal retina (Gomes ; He ).

Taking this information into consideration, we propose a model for group 1 retinoblastomas in which these tumours arise from a RPC in which the balance between P/P, P/D and D/D divisions is disrupted to favour continued proliferation, but with individual tumour cells showing stochastic regulation of lineage progression and variable patterns of gene expression. The higher level of apoptosis observed in group 1 retinoblastomas, compared with group 2, may reflect a proportion of tumour cells that undergo terminal D/D divisions. This model is consistent with published expression data on single retinoblastoma cells (McEvoy ). Examination of this data indicates that single cells obtained following repeated ( × 3) passage of a human retinoblastoma in mice do not show selection for a single retinal lineage but instead display a high level of inter-cell variability in the pattern of transcription factor expression (Supplementary Figure S5). These results support the idea that stochastic regulation of cell fate is an inherent characteristic of these retinoblastomas, as well as of normal retina.

Group 2 retinoblastomas in contrast are characterised by a stable pattern of cone photoreceptor gene expression, including not only genes encoding cone-enriched transcription factors but also multiple genes required for photoreceptor structure and function (Supplementary Table S2). These retinoblastomas might arise in a RPC in which the stochastic regulation of cell fate is strongly biased towards a cone cell lineage. Alternatively, they may originate in an RXRG-expressing precursor cell, primed for cone differentiation, which has failed to undergo cell cycle exit. Continued proliferation may be a consequence of the combined effects of pRB loss and also upregulation of genes responsible for cell cycle progression, for example, cyclin E2 and the mitotic cyclin-dependent kinase, CDK1, which showed 1.9- to 2.1-fold increased expression in group 2 relative to group 1 retinoblastoma and even larger increases relative to fetal and adult retina (results not shown).

The molecular aetiology of group 1 retinoblastomas

The restricted occurrence of chromosome 1q and 6p gain and 16q loss in group 1 retinoblastomas suggests that genes on these chromosomes are relevant to the aetiology of this tumour group. In agreement with previous reports, this study identified several functionally related genes with roles in spindle microtubule organisation and kinetochore/centrosome-mediated mitotic functions, for example, KIF14, CENPF, ASPM, NEK2 and NUF2 on chromosome 1 and KIFC1 on chromosome 6p (Grassmann ; Gratias ; Orlic ; Bowles ; Ganguly and Shields, 2010). These genes have particular relevance for group 1 retinoblastomas, as orientation of the mitotic spindle has been shown to influence the distribution of cell-fate determinants to one or both daughter cells, thereby producing symmetric (P/P, D/D) or asymmetric cell divisions (P/D or D/D with different cell types). In the retina, RPCs with horizontally oriented spindles (relative to the retinal pigment epithelium) generate identical daughter cells, while vertically oriented spindles tend to generate two different daughters, and it has been suggested that even slight changes in spindle orientation may influence the outcome of cell division (Cayouette ). Proteins such as ASPM have been shown to be critical for the maintenance of symmetric divisions in the developing CNS (Nigg and Stearns, 2011), and it is tempting to speculate that increased expression of ASPM might have a similar function in retinoblastoma, resulting in increased production of RPCs.

Genes on chromosomes 1q, 6p and 16q may also contribute to tumorigenesis through more direct effects on cell cycle regulation. An important consequence of loss of pRB function is aberrant regulation of the G1–S transition and these genes may further potentiate this effect through over-riding downstream cell cycle checkpoints. Of particular interest is the upregulation of the chromosome 1q gene CKS1B, which has an essential role in the SKP2-mediated ubiquitination of a range of cell cycle inhibitory proteins including p27/CDKN1B and also p130/RBL2 (Tedesco ; Krishnan ; Liberal ). CKS1B maps to 1q21, a common region of gain in cancer, and its amplification/overexpression has been reported in many tumours (Wang ; Chen ; Krishnan ). Also of interest is the downregulation of the 16q gene USP10, which de-ubiquinates p53 and is an essential regulator of its stability (Yuan ). It has been suggested that loss of USP10 provides an alternative mechanism to antagonise p53 function in tumours such as renal cell carcinoma, which (like retinoblastoma) do not have a high frequency of p53 mutations (Yuan ). USP10 was also identified by Gratias in a search for genes mapping to 16q24, which were differentially expressed in retinoblastomas with 16q LOH.

Although investigation of the molecular consequences of chromosome alterations in group 1 retinoblastomas provides clues about pathways which are aberrantly regulated in tumorigenesis, it is also important to identify the primary mitogenic signal(s). Potential candidates include growth factors, neurotrophins and neurotransmitters, all of which have been implicated in the regulation of cell proliferation in the developing retina (Martins and Pearson, 2008). Our GSEA analysis points to a role for neurotransmitter signalling, especially through amine and peptide ligand-binding receptors. This is a significant result as it provides further support for the novel finding of McEvoy ) that among paediatric solid tumour cell lines, the retinoblastoma cell lines Y79 and WERI1 are uniquely sensitive to broad-acting monoamine receptor inhibitors, and that blockade specifically of amine receptors (but not other major neurotransmitter signalling pathways) reduced retinoblastoma growth and survival. Inhibition of neurotransmitter signalling may therefore have important therapeutic potential in retinoblastoma.

The molecular aetiology of group 2 retinoblastomas

Group 2 retinoblastomas were observed to have a gene expression signature strongly indicative of differentiation along a cone photoreceptor lineage. This is consistent with the observation that human differentiated cone photoreceptors (but not early cone precursors) stain strongly for pRB, suggesting that pRB does not mediate the initial proliferative arrest in these cells (Lee ), but may be important for maintaining cell cycle exit and differentiation. A consequence of loss of RB function may be to permit cell cycle re-entry of cells, which otherwise are competent to form mature photoreceptors and have all of the components required for the visual cycle.

Our results suggest that there are two main features of cone photoreceptors which have relevance for tumorigenesis in group 2 retinoblastomas. The first is the extremely high metabolism of this cell type – photoreceptors are among the most metabolically active cells in the body (Reidel ). Of particular note in group 2 retinoblastomas was the upregulation of multiple components of the TCA and oxidative phosphorylation (OXPHOS) pathways. Although tumour cells are frequently reported to use aerobic glycolysis in preference to OXPHOS (Warburg effect), there is increasing evidence that OXPHOS and glycolysis cooperate to sustain energy needs in tumorigenesis, with dynamic regulation reflecting the tumour microenvironment (Jose ). Our results, and those of Chen et al (2007), suggest that in environments such as the retina and the brain where energy demand is high, tumour cells use enhanced mitochondrial respiratory pathways for energy production. The very high metabolic activity of retinal photoreceptor cells is also paralleled by elevated transcriptional and translational activity, and it is of interest that mutations of several pre-mRNA processing genes (e.g., PRPF3, PRPF31, PRPC8) are associated with retinal disease phenotypes (Cao ; Tanackovic ). Elevated mRNA processing, observed in group 2 retinoblastomas, may impact on tumorigenesis both through an overall promotion of cell growth and division, but potentially also through a generalised loss of splicing fidelity which is a hallmark of cancer (Ghigna ). Overall these results suggest that cellular and metabolic adaptations for efficient photoreceptor function in the normal retina may also contribute to a growth advantage for tumour cells in this microenvironment, and that this may account for the particular sensitivity of cone photoreceptor cells to oncogenic transformation.

Conclusions

In summary, the results of this study suggest that the molecular aetiology of retinoblastoma reflects aberrant regulation of developmental processes required for the highly ordered production of the different retinal cell types from a single progenitor cell. Although we cannot exclude the possibility that retinoblastomas represent a broad spectrum of developmental phenotypes, nevertheless our results suggest that most tumours arise from either a RPC or a committed cone precursor. Our findings are important for the development of novel targeted therapeutics for retinoblastoma and, in addition, strongly suggest that the clinical significance of retinoblastoma molecular heterogeneity merits further investgation.