PRACE ORYGINALNE/ORIGINAL PAPERS

Unsupervised analysis of follicular thyroid tumours transcriptome by oligonucleotide microarray gene expression profiling

Analiza nienadzorowana transkryptomu nowotworów pęcherzykowych tarczycy za pomocą profilowania ekspresji genów metodą mikromacierzy oligonukleotydowych

Bartosz Wojtaś^1*, Aleksandra Pfeifer¹, ^5*, Michał Jarząb^2*, Agnieszka Czarniecka³, Jolanta Krajewska¹, Michał Świerniak¹, Tomasz Stokowy¹,⁵, Dagmara Rusinek¹, Monika Kowal¹, Jadwiga Żebracka-Gala¹, Tomasz Tyszkiewicz¹, Małgorzata Oczko-Wojciechowska¹, Ewa Stobiecka⁴, Dariusz Lange⁴, Ralf Paschke⁶, Barbara Jarząb¹

¹Department of Nuclear Medicine and Endocrine Oncology, Maria Skłodowska-Curie Memorial Cancer Centre and Institute of Oncology, Gliwice Branch, Poland

²Department of Radiotherapy and Chemotherapy, Maria Skłodowska-Curie Memorial Cancer Centre and Institute of Oncology, Gliwice Branch, Poland

³Department of Oncologic Surgery, Maria Skłodowska-Curie Memorial Cancer Centre and Institute of Oncology, Gliwice Branch, Poland

⁴Department of Tumour Pathology, Maria Skłodowska-Curie Memorial Cancer Centre and Institute of Oncology, Gliwice Branch, Poland

⁵Faculty Of Automatic Control, Electronics And Computer Science, Silesian University of Technology, Gliwice, Poland

⁶Division of Endocrinology and Nephrology, University of Leipzig, Germany

^*All three authors (BW, AP, MJ) contributed equally

Barbara Jarząb, M. D., Ph. D. Department of Nuclear Medicine and Endocrine Oncology, Maria Skłodowska-Curie Memorial Cancer Centre and Institute of Oncology, Gliwice Branch, Wybrzeże Armii Krajowej 15, 44–101 Gliwice, Poland, tel.: +48 32 278 9339, fax: +48 32 278 9310, e-mail bjarzab@io.gliwice.pl

Introduction

Follicular thyroid cancer is a malignant tumour believed to arise from follicular cells in the thyroid gland. Being relatively rare compared to the papillary histotype of thyroid carcinoma, follicular cancer has not yet been thoroughly investigated in regard to the molecular mechanisms driving the malignant potential [1–4].

It is of special interest, because follicular thyroid carcinoma has a benign counterpart follicular thyroid adenoma. What is relatively unique, although resulting from the high degree of differentiation in this type of cancer and occurring also in other endocrine malignancies, is that there are no cell-specific criteria differentiating the areas with malignant potential from the fully benign tumours [5, 6]. The only defining feature is the presence of invasion of tumour capsule and vessels, if it is prominent enough to be spotted by a pathologist [7, 8]. The degree of invasion defines the subtypes of disease (minimally invasive, classical, and widely invasive) and certainly influences the reliability of the assessment by a pathologist. It is possible that some thyroid adenomas are in fact follicular carcinomas with no invasion revealed by pathology examination [4].

The difficulties described above give rise to the question whether the ability of follicular cells to invade is an inherent part of their genomic programme. If so, detection of that profile in cells could lead to genomic tests which on the basis of gene expression (and possibly in future also gene mutations and other aberrations) provide the information whether the tumour is benign or malignant not relying on the presence or absence of apparent invasion areas [9]. Unfortunately, comparisons on the global scale, i.e. between the full transcriptomes, almost always provides certain differences and their validation requires large independent samples [4, 9]. However, when the genomic distance between groups is large, the positive validation in further prospective studies is much easier and often leads to clinically significant conclusions.

One method to verify genomic distance between classes of samples is so called unsupervised analysis. In our bioinformatic study, we discriminate specimens on the basis of differences which are the major source of variability in the analysed dataset and are evident without previous assumptions regarding the sample class membership. In our study, we undertake an unsupervised analysis of a set of follicular thyroid tumours (both adenomas and carcinomas) to verify whether tumour malignancy is the factor influencing major sources of variability in our dataset.

Material and methods

We analysed a dataset consisting of 82 samples, among them thyroid carcinomas (FTC), thyroid adenomas (FTA) and grossly normal thyroid tissue (N) from the lobe of thyroid opposite to the tumours. Samples were fresh-frozen collected after surgery in the following Polish and German centres: MSC Memorial Cancer Centre and Institute of Oncology in Gliwice, University of Leipzig, University of Halle, and Mainz University Hospital. The study was approved by the local ethics committees and informed consent was obtained from all patients. The core set of samples (after exclusion of technical outliers in genomic study) consisted of 52 tumours (27 FTC, 25 FTA), of which 23% were collected from men. When histological re-assessment of all samples was carried out, two specimens were re-classified as a poorly differentiated thyroid carcinoma and one as a follicular variant of papillary thyroid cancer. Nine samples, among them seven carcinomas and two adenomas, were considered as oxyphilic tumours.

Total RNA was extracted by RNEasy columns (QIAGEN). Samples were analysed by oligonucleotide microarray HG-U133 Plus 2.0 according to the manufacturer’s protocol (Affymetrix). The detailed methods are given in a parallel manuscript describing supervised analysis of the dataset (Pfeifer et al., submitted, [10]). The dataset was pre-processed with the RMA and GC-RMA methods, routine R/Bioconductor packages were applied. Principal Component Analysis (PCA) was applied as a main method of unsupervised analysis. The dataset for this analysis was centered, no scaling was applied. We also used our own procedure of unsupervised gene selection based on Singular Value Decomposition of the dataset [11]. The results of unsupervised analysis are so called Principal Components which could be visualised on X-Y plots, showing two principal components on each plot. Samples might be coloured according to the class membership, but this step is done after the procedure was carried out. The lists of genes were analysed in regard of their Gene Ontology content and hierarchical clustering based on the lists was presented. The MSigDB database (Molecular Signatures Database, Broad Institute, MIT) was also used to analyse the biological meaning of selected gene lists. For each mode, an analysis of significant gene sets, as defined by MSigDB repository and with the use of GSEA algorithm was carried out. We also performed the analysis of an independent collection of normal human tissues, as provided by GNF Institute (Genomics Institute of Novartis Research Foundation). This was done by clustering of a selected subset of GNF data (genes for each mode), and in parallel displaying the same clustering for our samples.

Results

Unsupervised analysis of the full dataset

The first step was to apply Principal Component Analysis to all the samples obtained in the study. Principal Components #1 and #2 are shown on Fig. 1A. It is evident that certain samples (mainly carcinomas, but also three adenomas) are clear outliers in gene expression profile (they have clearly different profile to the majority of samples in the analysed group). Some of these effects might be influenced by the biological differences in these samples, but probably they are partially related to the technical issues (e.g. the extent of RNA degradation). When the same analysis was repeated after exclusion of 11 outlier samples (Fig. 1B), the clear pattern of differences between “normal” thyroid tissue (sample from the lobe contralateral to the tumour, depicted in red) and thyroid adenomas as well as carcinomas (depicted in green and blue) was seen.

A

B

Figure 1. Unsupervised analysis of the whole dataset (A) and after exclusion of 11 outlier samples (B). First two principal components are presented (PC#1 on X axis and PC#2 on Y axis). Normal thyroid samples (red dots), follicular thyroid adenomas (FTA, green) and follicular carcinomas (FTC, blue) appear in three overlapping clusters

Rycina 1. Nienadzorowana analiza całego zbioru danych (A) oraz zbioru danych po wyłączeniu 11 próbek odstających (B). Pierwsze dwie główne składowe przedstawiono na wykresie (składową 1 — PC#1 na osi X, składową 2 — PC#2 na osi Y). Próbki pochodzące z utkania tarczycy niezmienionej nowotworowo zaznaczono na czerwono, gruczolaki pęcherzykowe tarczycy na zielono, a raki pęcherzykowe tarczycy na niebiesko. Wyraźne jest znaczne nakładanie się tych trzech grup badanych

First and third Principal Components discriminate between carcinomas and adenomas

In the next step, two components which best discriminated between FTC and FTA, i.e. PC#1 and PC#3, were used together to provide unsupervised insight into differences between both groups (Fig. 2). Similarly to the previous observation, there was a trend towards differences in gene expression profile of follicular carcinomas (high values of PC#1, low of PC#3) and follicular adenomas (low values of PC#1, high of PC#3). When the artificial boundary between both clusters was arbitrarily drawn (black line), still the significant number of samples were misclassified. In the case of the adenoma cluster, at least four carcinomas were dispersed within, with one additional lying close to the group boundary. Within the cluster of carcinomas, at least three adenomas were indistinguishable in an unsupervised manner, with additional three within the ‘grey zone’. With the assumption of the boundary set by us, the unsupervised classification was concordant with the class membership in 76% of samples in the ‘left’ cluster dominated by carcinomas and in 79% of ‘right’ cluster samples, dominated by adenomas. These relatively high estimates confirm the significant genomic differences between two classes, simultaneously pointing out the outlier samples.

Figure 2. Unsupervised analysis of adenoma (FTA, red) and carcinoma samples (FTC, blue) by PC#1 and PC#3, outlier samples were excluded. Potential boundary between two clusters is marked by black line, with the number of outliers provided

Rycina 2. Analiza nienadzorowana zbioru gruczolaków pęcherzykowych (FTA, zaznaczone na czerwono) oraz raków pęcherzykowych (FTC, zaznaczone na niebiesko). Wykorzystano dwie główne składowe, składową pierwszą (PC#1) i trzecią (PC#3), wyłączono próbki znacznie odstające. Zaznaczono możliwą granicę między dwoma klastrami (czarna ciągła linia), przedstawiono także liczbę próbek nie poddających się takiej próbie klasyfikacji

Unsupervised analysis of all adenomas and carcinomas in the context of tumour histology

In the third step, we used the full set of 52 carcinomas and adenomas, with no exclusions (Fig. 3). Again, we might notice the difference between two populations when PC#1 and PC#3 are considered, although less visible than in the previous figure. When the information about histological parameters of samples was used to colour the points, it became evident that some of the outlier samples showed oncocytic features of the tumour, while the majority of non-oncocytic tumours cluster were within the main group of tumours. Also, two poorly differentiated tumours showed certain features of similarity to the group of outliers, although not evident enough to be identified in an unsupervised manner. Of interest, the one follicular variant of PTC which appeared in the analysed group did not differ strongly from the remaining follicular tumours in the context of PC#1 and PC3# variability.

A

B

Figure 3. Analysis of the core set of samples (without exclusion of outliers) by Principal Component Analysis. PC#1 and PC#3 are shown. A. It is apparent that neither PC#1 nor PC#3 provides the differentiation between carcinomas and adenomas; it is necessary to use the information from both principal components. B. When the same plot was coloured according to the oxyphilic features of tumours (red), some of these tumours appear to be the outliers with very low values of PC#1. Two poorly differentiated carcinomas (black) are also showing the relatively low values of PC#1. One sample inadvertently used in the study (follicular variant of papillary thyroid carcinoma, fvPTC, green) is not significantly outlying from the remaining follicular tumours

Rycina 3. Analiza głównego zbioru próbek (tylko raki i gruczolaki pęcherzykowe, ale bez wyłączania próbek odstających) za pomocą analizy głównych składowych. Przedstawiono dwie główne składowe, składową pierwszą (PC#1) i trzecią (PC#3). A. Wyraźnie widać że ani składowa pierwsza ani trzecia samodzielnie nie pozwalają na rozróżnienie obu grup, trzeba wykorzystać obie składowe równocześnie. B. Kiedy wykorzystano informacje dot. badań histopatologicznych próbki i zaznaczono te guzy w których stwierdzano cechy oksyfilności (czerwone); część z tych guzów ma charakter silnie odstający, z niskimi wartościami PC#1. Dwa niskozróżnicowane raki tarczycy (czarne) również wykazują względnie niskie wartości PC#1 w porównaniu z resztą próbek. Jedna z trzech próbek reklasyfikowanych (rak brodawkowaty wariant pęcherzykowy) wykazuje profil względnie zbliżony do pozostałych guzów

The biological meaning of observed differences

We carried out an extensive analysis of the biological character of genes correlated to the first six Principal Components. When genes, correlated to the first PC, were used to cluster carcinoma and adenoma samples, they appeared in two clusters: one with larger homogenous expression of subset of genes, relatively enriched in adenomas, and the other with various degrees of down-regulation of these genes and heterogeneous up-regulation of smaller cluster of transcripts, see Fig. 4A; the second cluster shows the relative enrichment of carcinomas. One may note that the clustering by the major source of variability provided by the first Principal Component is related neither to the tumour batch nor to the hybridisation batch, thus, the analysed variability does not arise from the technical factors (Fig. 4A).

When the analysed clustering was shown parallel to the clustering by the same list of genes (presented in the same order) on the collection of healthy human tissues (made available by Genomics Institute of Novartis Research Foundation, Fig. 4B) it was revealed that genes highly up-regulated, mainly in adenomas, were the thyroid-specific transcripts, some of them up-regulated also in liver (deiodinases). In the opposite, the second cluster of genes, with heterogeneous expression in the dendrogram branch enriched in carcinomas, contained immunity-related transcripts like HLA-DR, DQ genes and may be a genomic indicator of tissue infiltration by lymphocytes and other immune system cells.

In a similar manner, we analysed genomic context of the next Principal Components. Immunity-related transcripts were identified as those influencing the first, third and sixth principal components. First and third principal components were similar in their genetic content, as they contained the immunityrelated genes and provided the best discrimination between carcinomas and adenomas; third PC genes were not as apparently thyroid-related as PC#1 genes. Within the PC#4 we identified a large number of genes related to mitochondrial function; we consider this component as directly related to the tumour oxyphilicity features [12]. Within the fifth principal component, there were a large number of genes related to cell cycle and proliferation. A full analysis of the data is provided as a web appendix to this article, available online at www.zmnieo.home.pl.

A

B

Figure 4. Clustering of samples and genes by the gene list of transcripts correlated to the first Principal Component (PC#1). A. Clustering of the set of follicular carcinomas (FTC, blue) and adenomas (FTA, green). Sample batch and hybridisation batch are provided by colour to demonstrate that these two are not behind the observed variability. Left cluster of samples is relatively enriched in adenomas, while the right branch of dendrograms contains a large proportion of carcinomas

Rycina 4. Klasteryzacja próbek i genów za pomocą listy transkryptów skorelowanych z pierwszą główną składową (PC#1). A. Klasteryzacja zbioru raków (FTC, niebieskie) i gruczolaków pęcherzykowych (FTA, zielone). Seria próbki oraz hybrydyzacji mikromacierzowej są również zaznaczone kolorem, nie są one ewidentnie związane z głównymi źródłami zmienności. Lewy klaster dendrogramu zawiera nieco więcej gruczolaków, podczas gdy prawa gałąź grupuje więcej raków pęcherzykowych

Discussion

In our unsupervised analysis of follicular thyroid tumours, we identified the following sources of variability: oncocytic features [12], cell cycle and proliferation genes [13], immunity-related transcripts and expression of thyroid-specific genes [14].

Most of these have been shown previously, except immunity-related transcripts. Putting these factors in order of importance is not a trivial task, although it seems evident that thyroid-specific genes are the largest and strongest factor appearing in our analysis. On the other hand, immune response seems to be the most multifactorial and complicated variable, influencing to a lesser extent at least three principal components of the analysis.

In our previous analysis of papillary thyroid cancer (PTC) and contralateral normal thyroid tissues, we applied unsupervised analysis to seek for the potential sources of variability [15]. We found that the difference between PTC and normal thyroid was prominent and could be seen as the major source of variance. We also noted that immunity-related transcripts were the other strong factor influencing the overall profile of gene expression. We hypothesised that the variable expression of these transcripts was both endogenous (coming from expression in thyreocytes) and exogenous (transcripts from the immune cells infiltrating the tumour).

In the study we found that at least in certain follicular tumours it is possible to apply unsupervised analysis and discriminate between FTC and FTA. Approx. 80% of tumours could easily be discriminated between clusters, while the remaining samples lie in between and could not be distinguished. It is certainly difficult to judge whether histopathological diagnosis of the analysed samples is accurate, and whether some of the follicular adenomas are in fact tumours with gene expression profile of cancer, but do not manifest the invasive properties. Also the opposite misclassification (falsely positive diagnosis of cancer) could influence results. Our previous studies on the accuracy of histopathological diagnosis of thyroid cancer [16] indicate that the misclassification rate is substantial, even among experienced histopathologists.

There have been a number of genomic studies dealing with the classification of follicular thyroid tumours by genomic tests [17–22]. The recent paper by Alexander et al. provides a clinically useful tool to discriminate gene expression patterns of thyrocytes even in fine-needle aspiration biopsy [23–24]. However, all the studies carried out until now do not provide the ideal classification accuracy for follicular tumours, and they claim a misclassification rate of between 5% and 30%. Thus, we undertook our study to verify the genomic distance between FTC and FTA and test whether it is large enough to provide good results of molecular classifiers. Our previous analyses [25] showed that for papillary thyroid cancer the genomic distance is high, and thus the misclassification rate might be very low. For follicular neoplasms analysed in the present study, we observed significant overlap between the two datasets and we estimate that the misclassification rate must be between 10% and 20% in an unsupervised approach.

Conclusions

In conclusion, the difference between follicular thyroid carcinoma and thyroid adenoma is moderate. Approx. 80% of samples could be discriminated in unsupervised analysis, while the remaining samples constitute the grey zone and could not be discriminated in that way. At least two principal components are necessary to provide such discrimination, confirming the multifactorial mechanisms driving the invasive potential of follicular carcinomas. Partially the expression of thyroid-specific genes and immunity-related transcripts improves the unsupervised classification of follicular tumours. The influence of oxyphilic character of some thyroid tumours might also be a strong factor appearing as a strong source of gene expression variability.