Clustering is een machine learning-techniek uit de categorie unsupervised learning. Hoe gaat het in z’n werk? Een voorbeeld uit de genetica, waarbij celtypes van bloedcellen achterhaald worden op basis van de genexpressie van de cellen, laat de stappen goed zien.

In een vorige baan heb ik gewerkt met single-cell RNA sequencing-data, ofwel genexpressiedata van individuele cellen. Dit is een vrij nieuwe techniek: waar eerder genexpressie werd gemeten van meerdere cellen tegelijk, kan dit nu voor elke cel apart gedaan worden. Dit kan interessante inzichten geven in de heterogeniteit tussen cellen, maar informatie over de celtypes gaat verloren in het proces.

Single-cell RNA sequencing

Een gen in het DNA komt tot uiting als het wordt omgezet in RNA, dat vervolgens als blauwdruk dient voor de bouw van een eiwit. Door naar het RNA in een cel te kijken, weet je dus welke genen tot expressie komen in die cel.

Bij single-cell RNA sequencing wordt een groep cellen zo bewerkt dat elke cel een unieke code krijgt die aan al het RNA binnen de cel wordt geplakt. Daarna wordt het RNA van alle cellen uitgelezen (gesequenced) door een machine. Met behulp van de unieke celcodes kan het uitgelezen RNA van elke cel achteraf weer bij elkaar gezocht worden. Het resultaat van dit proces is een tabel met de cellen als rijen en de genen als kolommen, die de gemeten genexpressie per cel weergeeft. De cellen zelf worden vernietigd in dit proces: als er verschillende celtypes tegelijk gesequenced worden, is niet bekend welke cel van welk type was.

De computer kan na het sequencen op zoek naar structuur in de tabel: zijn er groepen cellen te vinden die qua genexpressie sterk op elkaar lijken? Dit proces wordt clustering genoemd. Bij de genexpressiedata van verschillende soorten bloedcellen blijkt clustering geschikt om de celtypes terug te vinden1.

De data

Single-cell RNA sequencing is onder andere toegepast op PBMCs (peripheral blood mononuclear cells). Dit is een verzamelnaam voor verschillende types bloedcellen met een ronde kern, bestaande uit lymphocyten (T-cellen, B-cellen en NK-cellen) en monocyten. Een sequencing-dataset van 2.700 PBMC-cellen van een gezond persoon is online beschikbaar op de website van 10x Genomics, een fabrikant van single-cell RNA sequencing machines2. Online zijn ook uitgebreide tutorials in R en Python te vinden om met clustering de celtypes van deze dataset te achterhalen3.

Clustering stap voor stap

Stap 1: data voorbereiden

In de eerste stap van een clusteringproject gaat het om de dataset verkennen en de kwaliteit in kaart te brengen. Ook is vaak een vorm van datanormalisatie nodig om de getallen in de dataset te standaardiseren.

De dataset bestaat uit 2.700 cellen en 32.738 genen. Bij het sequencen gaat data verloren: niet al het RNA wordt goed opgepikt. Er zullen daardoor cellen zijn waarvan maar heel weinig RNA is uitgelezen. Deze cellen kunnen het beste uit de dataset verwijderd worden. Ook genen die maar in heel weinig cellen geteld zijn, kunnen beter buiten beschouwing worden gelaten. Na nog enkele andere kwaliteitschecks blijven er 2.638 cellen en 13.714 genen over.

Niet elk gen is even interessant: het doel is om celtypes te onderscheiden door te kijken naar verschillen in genexpressie. Genen die in alle cellen ongeveer evenveel tot expressie komen, zijn daarom veel minder interessant dan genen die in sommige cellen veel en in andere cellen juist weinig tot expressie komen. Door te filteren op dit soort variabele genen, blijven er nog 1.838 genen over. De dataset is daarmee dus een stuk kleiner geworden. Dat is goed nieuws voor de vervolgstappen, maar het is nog niet genoeg.

Stap 2: de dimensionaliteit verlagen

Met nog 1.838 genen in de dataset, levert elke cel een 1.838-dimensionaal datapunt op. Veel algoritmes hebben moeite met dimensies die zo hoog zijn, daarom moet de dimensionaliteit eerst verder gereduceerd worden.

Een manier om dat te doen is met principale-componentenanalyse (PCA). PCA kan de belangrijkste variatie in een dataset samenvatten in zogeheten principale componenten. Van elke van deze componenten is uit te rekenen hoeveel variatie ze beschrijven. Het blijkt dat een handjevol van deze principale componenten al gauw het belangrijke deel van de variatie in een dataset kan vatten. Dat is ook het geval bij deze dataset, waar de eerste 7 principale componenten relatief veel variatie verklaren:

De verklaarde variantie per principale component van de bloedcellen daalt snel na de 7e component

Na de zevende component knikt de lijn in de grafiek, en vanaf daar voegt elke volgende component nog maar weinig extra toe. De genexpressie van elke cel kan samengevat worden in termen van deze eerste 7 principale componenten, waardoor voor elke cel slechts een 7-dimensionaal datapunt overblijft.

Stap 3: de dataset visualiseren

De vorige stap leverde een flinke verlaging op van de dimensionaliteit, en de data is nu geschikt om te gaan clusteren. Maar een 7-dimensionale dataset heeft nog wel het nadeel dat je hem niet in z’n geheel kan visualiseren. Daarom zijn er technieken om de structuur van de data in deze hogere dimensie te projecteren naar twee dimensies.

Voorbeelden hiervan zijn t-SNE en UMAP. Het zijn beide stochastische technieken, wat wil zeggen dat het resultaat kan verschillen als het algoritme meerdere keren gedraaid wordt op dezelfde dataset. Ook geven de geprojecteerde waardes op de x– en y-as geen goede weergave van de daadwerkelijke afstand tussen punten, daarom worden de waardes op de assen vaak weggelaten.

Het resultaat van t-SNE toepassen op de bloedcellen geeft de volgende figuur. Elk punt stelt een cel voor, en cellen die op elkaar lijken worden dicht bij elkaar geprojecteerd. Meerdere groepen worden nu zichtbaar:

de tSNE-plot van de bloedcellen laat al enkele groepen zien

Stap 4: clustering

Nu kan er geclusterd worden op de principale componenten. Er bestaan veel verschillende clusteringalgoritmes en in Scikit-learn zijn er veel geïmplementeerd4. Een van de eenvoudigste is K-means. Dit algoritme begint met een aantal willekeurige centra, en wijst vervolgens elk datapunt toe aan het dichtstbijzijnde centrum. Vervolgens worden de centra opnieuw berekend als het gemiddelde van alle toegewezen datapunten en begint het toewijzingsproces opnieuw. Het proces stopt als de toewijzing van punten aan de centra niet meer verandert.

Bij K-means moet je zelf kiezen hoeveel clusters je wilt hebben. De keuze voor het aantal clusters kun je op verschillende manieren maken. Een daarvan is door de silhouette score te berekenen voor verschillende aantallen clusters. De silhouette score combineert een maat voor hoe dicht elk datapunt bij andere punten in zijn eigen cluster ligt met hoe ver het verwijderd is van punten in het dichtstbijzijnde andere cluster. De score zit tussen 0 en 1, en hoe hoger de score hoe beter de clustering. Het blijkt dat in dit geval de silhouette score bij 8 clusters het hoogste is:

Silhouette score voor verschillende aantallen clusters laat zien dat clustering met 8 clusters de hoogste score geeft

Als de tSNE-afbeelding wordt gekleurd op basis van deze 8 clusters, ziet dat er als volgt uit:

Clustering toegepast op de bloedcellen laat binnen de groepen van de tSNE plot nog weer subgroepen zien

Dat begint ergens op te lijken, maar de cruciale stap moet nog komen. Welke cluster komt nu overeen met welk celtype?

Stap 5: betekenis geven aan de clusters

Nu er clusters gevonden zijn, moet nog gekeken worden of de clusters wel betekenisvol zijn. Wat voor groepen zijn er te verwachten in de data? Welke kenmerken onderscheiden die groepen?

De cellen die we hier bekijken, zijn PBMCs. Over deze celtypes is veel bekend. Zo weet men dat het MS4A1-gen vooral in B-cellen tot uiting komt, terwijl het gen GNLY vooral tot expressie komt in NK-cellen. Op dezelfde manier heeft elk celtype kenmerkende genen.

Hieronder is de tSNE-plot gekleurd op basis van het MS4A1-gen (links) en GNLY (rechts). Hoe donkerder elke cel kleurt, hoe sterker het gen tot expressie komt in die cel:

De kenmerkende genen komen duidelijk in verschillende clusters tot expressie. Het lijkt er dus op dat de verschillende clusters inderdaad celtypes voorstellen. De eerste twee clusters kunnen nu geïdentificeerd worden met hun waarschijnlijke celtype:

Celtypes toekennen aan de clusters die volgden na de clustering

Door verder in te zoomen op alle kenmerkende genen en hun expressie binnen de verschillende clusters, kunnen andere celtypes geïdentificeerd worden. Er blijken verschillende soorten T-cellen en monocyten te kunnen worden onderscheiden (zie ook de eerdergenoemde tutorials).

Kwaliteit van de clusters

Hoewel het GNLY-gen in het cluster linksboven heel duidelijk naar voren kwam, zijn er ook enkele cellen uit andere clusters met een hoge expressie van GNLY. Klopt het dan wel dat dit geen NK-cellen zijn? Dat is het lastige bij clusteren. Waar hoort de grens tussen clusters precies te liggen? Verschillende algoritmes geven verschillende antwoorden. Hoe accuraat de uiteindelijke celtypebepaling door middel van clustering is, is op basis van alleen deze dataset niet te achterhalen.

Ook is het belangrijk om in het achterhoofd te houden dat de meeste algoritmes elk punt aan een cluster toekennen, ook als dat punt eigenlijk niet zo goed bij dat cluster past – het past alleen nog minder goed bij de andere clusters. En soms zijn er helemaal geen betekenisvolle clusters te vinden, of blijken de clusters niet betekenisvol te zijn voor het beoogde doel.

Tot slot

Om clustering toe te kunnen passen, moet de data dus eerst worden voorbereid. Vaak moet de dimensionaliteit van de dataset ook worden verlaagd: in dit voorbeeld ging de dataset van bijna 33.000 genen naar slechts 7 principale componenten. Daarmee kon al een redelijke clustering gevonden worden!

Er zijn verschillende clusteringsalgoritmes, die verschillende resultaten geven. Vervolgens moet je wel zelf betekenis geven aan de gevonden clusters, waar domeinkennis (in dit geval kennis over kenmerkende genen van celtypes) voor nodig is.

Al met al kan clustering een krachtige techniek zijn om structuur te vinden in een dataset als de labels – de antwoorden van wat je wilt voorspellen, in dit geval informatie over de celtypes – ontbreken. Als deze labels er wel zijn, kunnen supervised learning-algoritmes uitkomst bieden. Hier lees je meer over de verschillen tussen supervised en unsupervised learning.

Bronnen

  1. Satija, R., Farrell, J., Gennert, D. et al. Spatial reconstruction of single-cell gene expression data. Nat Biotechnol 33, 495–502 (2015). DOI ↩︎
  2. Dataset: http://cf.10xgenomics.com/samples/cell-exp/1.1.0/pbmc3k/pbmc3k_filtered_gene_bc_matrices.tar.gz ↩︎
  3. Guided clustering tutorial in R (website Satijalab) en PBMC-tutorial met scanpy (Python) ↩︎
  4. Documentatie verschillende clustering-algoritmes in Scikit-learn ↩︎

Een Python-notebook van bovenstaande stappen is te vinden op GitHub.