sejttípusok azonosítása egysejtű adatokból stabil klaszterezés segítségével

a javasolt módszer célja az egyes sejtek keverékében jelen lévő sejttípusok azonosítása. A módszer bemenete az egysejtű génexpressziós mátrix (Mgene), amelyben a sorok a géneket, az oszlopok pedig a sejteket képviselik. Az alábbiakban részletesebben ismertetjük a bemeneti adatokat és a javasolt keret különböző lépéseit. Az általános megközelítés ábrán látható. 1.

1.ábra

a javasolt módszer teljes munkafolyamata. Tekintettel az egysejtű génexpressziós mátrixra, az (A) modul kiküszöböli azokat a géneket, amelyek egyetlen sejtben sem expresszálódnak. A kapott mátrix felhasználásával a (B) modul kiszámítja a sejtek közötti euklideszi távolságot. Ennek a modulnak a kimenete egy távolságmátrix, amelyben a sorok és oszlopok a cellák(dcell). A (C) modul csökkenti a távolságmátrix dimenzióját a t-elosztott sztochasztikus szomszéd beágyazási (t-SNE) technikával. Ebben a modulban egy átlagos sziluett módszert alkalmaznak a klaszterek optimális számának kiválasztására k. végül a (D) modulban az alsó dimenziós távolságmátrixot és a (C) modulból kapott k klaszterek optimális számát használják bemeneti adatként a cellák legstabilabb csoportosulásának azonosításához. A 2. ábra a D modul részleteit mutatja.

adatforrás

a nyolc nyilvánosan elérhető scRNA-seq adatkészlet, valamint az elemzésünkben használt öt szimulációs adatkészlet szerepel a kiegészítő anyagokban. A nyolc valós adatkészlet közül három kivételével az összes (Klein51, Patel52, Treutlein53) arany standardnak tekinthető, mivel a cellák címkéi véglegesek. A Patel52-et és a Treutlein53-at Kiselev et al.28 mivel sejtcímkéiket a számítási módszerek és a szerzők által az alapul szolgáló biológiáról szerzett ismeretek alapján határozzák meg.

a feldolgozott adatokat a Hemberg lab weboldaláról szereztük be (https://hemberg-lab.github.io/scRNA.seq.datasets). Hemberg et al.54 használja a SingleCellExperiment Bioconductor S4 class55-öt az adatok tárolására, a scater package56-ot pedig minőségellenőrzési és ábrázolási célokra. A normalizált adatok letétbe, mint egy SingleCellExperiment objektum (.RData fájl) és a cellatípus információ az objektum “colData” nyílásának cell_type1 oszlopában érhető el. A sejtek génexpressziós értékei mátrixként vannak rendezve, amelyben a sorok sejtek, az oszlopok pedig a gének. Elemzésünk során eltávolítják azokat a géneket (jellemzőket), amelyek egyetlen sejtben sem expresszálódnak. Ebben az elemzésben egyetlen cellát sem szűrtünk.

Génszűrés

Az ábrán látható módon. 1A, eltávolítjuk azokat a géneket/transzkriptumokat, amelyek egyetlen sejtben sem expresszálódnak (az expressziós érték minden sejtben nulla). Az ilyen gének nem tudnak olyan hasznos információkat szolgáltatni, amelyek különbséget tudnak tenni a sejttípusok között57. A szűrési módszer egysejtű génexpressziós mátrixon (Mgene) történő végrehajtásának eredményét használjuk a javasolt keretrendszer második moduljának bemeneteként.

A sejtek közötti különbség mérése

a sejtek közötti távolságot az euklideszi metrika segítségével számítjuk ki (ábra. 1B). Ennek a lépésnek a kimenete a dcell dcell (dissimilarity) mátrix dcell (dcell) cella. A D dimenzióját a t-elosztott sztochasztikus szomszéd Beágyazással (t-SNE)34,58, a nemlineáris dimenzió redukciós/vizualizációs technikával (ábra. 1C). Akkor olvassa el a kimeneti D’cell×l, ahol 2 ≤ l ≤ sejt. Ebben a tanulmányban a méretek száma 2.

klaszterezés

a klaszterek optimális számának meghatározása

Ez a rész a javasolt módszer harmadik modulját ismerteti (ábra. 1C). Ebben az elemzésben a t-SNE-t ismételten (n = 50) alkalmazzák a dcell dcell távolságmátrixra (d’cell), hogy megkapjuk a D’ cell dimenziócsökkentett távolságmátrixot (d ‘cell). l.minden alkalommal a klaszterek optimális számát az átlagos sziluett módszer alapján számítjuk ki, a D’ dimenziócsökkentett távolságmátrix felhasználásával. Annak érdekében, hogy megtaláljuk a K klaszterek optimális számát, A k-azt jelenti, hogy a klaszterezést a D’ mátrixra alkalmazzuk egy tartományérték segítségével (alapértelmezett = 2:20), és kiválasztjuk azt a k-t, amely maximalizálja az átlagos sziluettmérést. Végül a kiválasztott számok átlagát k különböző ismétlések között (n = 50) (a legközelebbi egész számra kerekítve) a klaszterek végső optimális számának tekintjük.

a silhouette az adott csoportosítás minőségét annak alapján értékeli, hogy az adatpontok mennyire vannak csoportosítva. Az egyes adatpontokhoz sziluettmérés van hozzárendelve, amely azt mutatja, hogy egy adatpont milyen közel van a saját klaszteréhez a többi klaszterhez képest. Minden I adatpontra ezt a mértéket a következőképpen számítjuk ki:

$${\rm{s}}({\rm{i}})=\frac{b(i)-A(i)}{max\{a(i),b(i)\}}$$

ahol a(i) az I adatpont és az ugyanazon fürt összes többi adatpontja közötti átlagos távolság. b (i) az I legkisebb átlagos távolsága bármely más halmaz összes pontjától, amelynek i nem tagja. az s (i) -1-től 1-ig terjedő értékeket vesz fel, ahol a magas pozitív pontszám azt mutatja, hogy az adott adatpont jól csoportosul (közel a saját klaszterének más pontjaihoz, távol a többi klaszter pontjaitól). Ezzel szemben a magas negatív pontszám azt mutatja, hogy az adatpont rosszul csoportosul.

k-az újramintavételi módszeren alapuló klaszterezést jelenti

Ez a szakasz a javasolt módszer utolsó moduljának részleteit ismerteti. Amint az ábrán látható. 2. A D’ dimenzionalitás csökkentett távolság mátrix és az előző lépésből kiválasztott k klaszterek számának felhasználásával különböző klaszterezési megoldások generálásával azonosítjuk a legstabilabb klaszterezést (clusteringi (i)), és újramintázási módszerrel mérjük az egyes klaszterezési megoldások stabilitását. Az egyes klaszterekhez rendelt stabilitási intézkedés (klaszterezés) azt mutatja, hogy az adott klaszterhez tartozó k klaszterek milyen gyakran maradnak meg, ha a bemeneti adatokat (D’) többször újramintázzák. Az újramintázott adatkészleteket D’ – ből generálják úgy, hogy az adatpontok (cellák) 5% – át véletlenszerűen helyettesítik zajjal. Ezeket a zajos adatkészleteket ezután a K-means algoritmus bemeneteként használják. Ezért számos klaszterezés (clusteringi, j, J) keletkezik az újramintázott adatokból (a clusteringi újramintázott verziói).

2.ábra

a legstabilabb csoportosítás azonosítása. Ebben az elemzésben, figyelembe véve az alsó dimenziós d ‘ cell távolság mátrixot, az L-t és a K klaszterek optimális számát, kiszámítjuk n különböző klaszterek (klaszterezés1,…, klaszterezésn) a k – azt jelenti klaszterezési algoritmus. Ezután az egyes klaszterek stabilitását újramintázási megközelítés alapján értékelik (szürke doboz). Az egyes fürtözésekhez stabilitási pontszámot rendelnek annak alapján, hogy a klaszterek milyen gyakran kerülnek helyreállításra a bemeneti adatok zavarásakor (újramintázva). Végső megoldásként a maximális stabilitási pontszámmal rendelkező klaszter kerül kiválasztásra.

annak érdekében, hogy felmérjük az egyes C klaszterek stabilitását a klaszterezésbeni (eredeti klaszterezés), a C klasztert összehasonlítjuk a klaszterezés összes klaszterével, amelyet az újramintaadatokból (klaszterezési, j) kapunk a Jaccard távolság alapján. A Jaccard együttható59, a halmazok közötti hasonlóság mértéke, két klaszter hasonlóságának kiszámítására szolgál az alábbiak szerint:

$${\rm{J}} ({\rm{A}}, {\rm{B}})= \ frac {/a \ cap B/} {/A \ cup B/},\, A, B\subseteq x$$

ahol az A és B kifejezés két klaszter, amelyek az X = {x1,…, xN} néhány adatpontjából állnak.

Ha a Jaccard hasonlósága a C fürt (az eredeti fürtözési fürtözésből) és az újramintázott fürtözésben a leginkább hasonló fürt között egyenlő vagy nagyobb, mint 0,75, akkor a fürt stabilnak tekinthető (megőrzött). Így a stabilitás, az egyes klaszter clusteringi kell kiszámítani, mint a százalékos alkalommal, hogy a halmaz tartósítva (Jaccard együttható ≥ 0-ra.75) az egész m különböző resamplings.

ezután átlagoljuk a clusteringi-hez tartozó K klaszterek stabilitási mutatóit, és ezt tekintjük a clusteringi Általános stabilitási mutatójának. N különböző klaszterezési megoldás közül (clusteringi (i)) a maximális stabilitási mutatóval rendelkező klaszterezési megoldást választjuk végső klaszterezési megoldásnak.

a 3. ábra mutatja az egyes klaszterek stabilitási mértékének kiszámításához elvégzett resampling módszer részleteit. A K-átlag alkalmazásával kapott klasztereket az újramintázott adatkészleten csak a nem zajpontok alapján hasonlítják össze az eredeti bemeneti adatokból származó klaszterekkel (a zajadatpontok kizárásra kerülnek, ha két klasztert hasonlítanak össze a Jaccard hasonlósági mutató alapján.

3.ábra

az újramintavételi keretrendszer az egyes klaszterek stabilitási mértékének kiszámításához. A bemenet N adatpontot tartalmaz X = {x1,…, xN}, a klaszterek száma k, az újracsomagolások száma m, valamint a C klaszterezés, amelyet k-eszközök alkalmazásával kapunk X. Ez az elemzés generál m újramintavétel adatok az adatpontok 5%-ának véletlenszerű helyettesítésével a zajjal, és kiszámítja m újramintázott klaszterek alapján k-eszközök klaszterezés. A C-ben lévő c klasztereket összehasonlítjuk a resampling klaszterezés leginkább hasonló klaszterével, és kiszámítjuk a két klaszter közötti Jaccard-együtthatót, miközben a zajpontokat kizárjuk. Azoknak az időknek a százalékos aránya, amikor a Jaccard-együtthatók nagyobbak, mint 0.A 75. szám a C klaszter stabilitási mércéje. a C klaszterhez tartozó összes klaszter stabilitási mérésének átlagát számítjuk ki és tekintjük a C klaszterezés Általános stabilitási mércéjének.

validálási módszerek

13 különböző adatkészletet használunk, amelyekben a cellatípusok (címkék) ismertek. A referenciacímkék és az egyes fürtözési módszerekkel kapott kikövetkeztetett címkék hasonlóságának mérésére három különböző mutatót használunk: korrigált rand index (ARI), korrigált kölcsönös információ (AMI) és V-mérés az alábbiak szerint.

Korrigált rand index

a cellacímkék alapján a korrigált Rand indexet (Ari)47 használjuk a következtetett csoportosítás és a valódi csoportosítás közötti hasonlóság értékelésére. ARI tól 0, a rossz megfelelő (véletlenszerű klaszterezés), hogy 1 egy tökéletes megállapodás a valódi klaszterezés. N adatpontok halmaza esetén a kontingenciatábla a két klaszter közötti megosztott adatpontok száma alapján készül. Tegyük Fel, Hogy X = {X1, X2, …, XR} és Y = {Y1, Y2, …, YC} két különböző klasztert képvisel R, illetve C klaszterekkel. Az X és Y közötti átfedés egy kontingenciatáblázatban foglalható össze , ahol az MR = C=, ahol i = 1…R, j = 1…A C. Xi és Yj az X és Y csoportosulásokban egy klasztert jelöl, az i és j pedig a készenléti táblázat sorszámára, illetve oszlopszámára utal. Az ARI meghatározása a következő:

Korrigált kölcsönös információ

$$H(X)=\mathop{\sum }\limits_{i\mathrm{=1}}^{R}P(I)\,logP(i)$$
(2)

H(X) nem negatív, és csak akkor veszi a 0 értéket, ha nincs bizonytalanság az adatpont klasztertagságának meghatározásában (csak egy klaszter). A két X és Y csoportosítás közötti kölcsönös információ (Mi) kiszámítása a következőképpen történik:

$$MI(X,Y)=\mathop{\sum }\limits_{i\mathrm{=1}}^{R}\mathop{\sum }\limits_{j\mathrm{=1}}^{C}P(i,j)\,log\frac{P(i,j)}{P(i)P(j)}$$
(3)

ahol P(I, J) azt a valószínűséget jelöli,hogy egy adatpont mind az X X, mind az YJ klaszterhez tartozik y:

$$p(i, j)=\frac{|{x}_{i}\cap {y}_{j}|}{n}$$
(4)

mi egy nem negatív mennyiség, amelyet az entrópiák határolnak H(X) és h(y). Számszerűsíti a két csoportosulás által megosztott információkat, ezért csoportosítási hasonlósági intézkedésnek tekinthető. A kölcsönös információ korrigált mértéke a következőképpen van meghatározva:

$$AMI(X,Y)=\frac{MI(X,Y)-E\{Mi(X,Y)\}}{max\{H(X),H(Y)\}-e\{mi(X,Y)\}}$$
(5)

ahol a várható kölcsönös információ két véletlenszerű klaszter között:

ahol az AI és a bj a kontingencia táblázat részleges összegei: \({a}_{I}={\sum }_{J\mathrm{=1}}^{c}{n}_{IJ}\) és \({B}_{J}={\sum }_{I\mathrm{=1}}^{r}{n}_{IJ}\).

a korrigált kölcsönös információ (ami) értéke 1, ha a két klaszter azonos, és 0, ha a mi két partíció között megegyezik a véletlen miatt elvárt értékkel.

V-mérték

A V-mérték 50 a két mérték közötti harmonikus átlag: homogenitás és teljesség. A homogenitási feltétel akkor teljesül, ha egy klaszterezés csak azokat az adatpontokat rendeli hozzá egyetlen klaszterhez, amelyek egyetlen osztály (igaz fürt) tagjai. Így az egyes klasztereken belüli osztályeloszlást egyetlen osztályra kell torzítani (nulla entrópia). Annak meghatározásához, hogy egy adott csoportosítás milyen közel van ehhez az ideálhoz, az osztályeloszlás feltételes entrópiáját az azonosított csoportosulás alapján számítjuk ki H(C / K), ahol C = {C1, C2,…, Cl} osztályok halmaza, K pedig klaszterezés K = {K1, K2,…, Km}. Tökéletesen homogén esetben ez az érték 0. Ez az érték azonban az adatkészlet méretétől és az osztályméretek eloszlásától függ. Így ezt a feltételes entrópiát normalizálja az entrópia maximális csökkenése, amelyet a klaszterezési információ nyújthat, H (C). Ezért a homogenitást a következőképpen határozzuk meg:

$$h=\{\begin{array}{cc}1 & \text{if}\,H(C,K)=0\\ 1-\frac{H(C| K)}{H(C)} & \text{egyébként}\end{array}$$
(7)

a teljesség szimmetrikus a homogenitással50. A teljességi kritériumok teljesítése érdekében a fürtözésnek minden olyan adatpontot hozzá kell rendelnie, amely egyetlen osztály tagja, egyetlen fürthöz. A teljesség méréséhez felmérjük a klaszterfeladatok eloszlását az egyes osztályokon belül. Egy tökéletesen teljes fürtözési megoldásban ezen eloszlások mindegyike teljesen egyetlen klaszterre torzul.

a H homogenitás és a C teljesség alapján a V-mérték a homogenitás és a teljesség súlyozott harmonikus átlagaként kerül kiszámításra:

$${\rm{V}} \mbox{-} {\rm{m}}{\rm{e}}{\RM{S}}{\rm{u}}{\rm{r}}{\rm{e}}=\frac {(1+\beta)\AST h\AST C} {(\beta\AST H)+C}$$
(8)

ha a 6-nál nagyobb, akkor a teljességet erősebben súlyozzuk a számítás során. Ha az 1-nél kisebb, akkor a homogenitás erősebben súlyozódik. Mivel a homogenitás, a teljesség és a V-mérték számításai teljesen függetlenek az osztályok számától, a klaszterek számától, az adatkészlet méretétől és a klaszterezési algoritmustól, ezek a mérések bármilyen klaszterezési megoldás értékelésére alkalmazhatók.

Vélemény, hozzászólás?

Az e-mail-címet nem tesszük közzé.