Massivt parallell digital transkriptionsprofilering av enskilda celler

Höghastighetsavbildning av gelpärlor och celler i GEMs

Ett mikroskop (Nikon Ti-E, × 10-objektiv) och en höghastighetsvideokamera (Photron SA5, bildfrekvens=4 000 s-1) användes för att avbilda varje GEM när den skapades i det mikrofluidiska chipet. En skräddarsydd analysprogramvara användes för att räkna antalet genererade GEM och antalet pärlor i varje GEM, baserat på kantdetektering och kontrasten mellan pärlkanter och GEM-kanter och den intilliggande vätskan. Resultaten av analysen sammanfattas i figur 1c. För att uppskatta fördelningen av celler i GEMs användes manuell räkning för ∼28k bilder av en video på en delmängd GEMs. Resultaten visar på en ungefärlig följsamhet till en Poissonfördelning. Procentandelen av inkapslingar av flera celler var dock 16 % högre än det förväntade värdet, vilket möjligen kan bero på fel i delprovtagningen eller på interaktioner mellan celler och celler (några klumpar med två celler observerades under den manuella räkningen) (kompletterande fig. 1b).

Cellinjer och transplantationspatientprover

Jurkat (ATCC TIB-152), 293T (ATCC CRL-11268) och 3T3 (ATCC CRL-1658) celler förvärvades från ATCC och odlades enligt ATCC:s riktlinjer. Färska PBMC, frysta PBMC och BMMC köptes från ALLCELLS. Frysta PBMC från donator A tillverkades av färska PBMC från donator A genom att försiktigt blanda 1e6 celler i frysmedium (15 % dimetylsulfoxid (DMSO) i Iscoves modifierade Dulbecco-medium innehållande 20 % FBS) och kyla dem i CoolCell FTS30 (BioCision) i -80 °C i minst 4 timmar innan de överförs till flytande kväve för lagring i 3 veckor.

Institutional Review Board vid Fred Hutchinson Cancer Research Center godkände studien av transplantationsprover. De förfaranden som följdes var förenliga med Helsingforsdeklarationen från 1975 och den gemensamma regeln. Proverna togs efter att patienterna hade gett skriftligt informerat samtycke till molekylära analyser. Vi identifierade patienter med AML som genomgått allogen hematopoetisk stamcellstransplantation vid Fred Hutchinson Cancer Research Center. Diagnosen AML fastställdes enligt Världshälsoorganisationens reviderade kriterier33.

Benmärgsaspirat erhölls för klinisk standardtestning 20-30 dagar före transplantationen och seriellt efter transplantationen enligt behandlingsprotokollet. Alikvot av benmärgsaspirat bearbetades inom 2 timmar efter uttag. BMMCs isolerades med hjälp av centrifugering genom en Ficollgradient (Histopaque-1077; Sigma Life Science, St Louis, MO, USA). BMMC:s samlades upp från gränssnittet mellan serum och Ficoll med en engångspipett av Pasteur-typ och överfördes till ett 50 ml koniskt rör med 2 % patientserum i 1 × PBS. BMMC:erna räknades med en hemacytometer och livskraften bedömdes med hjälp av Trypanblått. BMMC:erna resuspenderades i frysmedium med 90 % FBS och 10 % DMSO och frystes med hjälp av en Thermo Scientific Nalgene Mr Frosty (Thermo Scientific) i en frys på -80 °C i 24 timmar innan de överfördes till flytande kväve för långtidsförvaring.

Beskattning av RNA-innehåll per cell

Mängden RNA per celltyp bestämdes genom att kvantifiera (Qubit; Invitrogen) RNA som extraherats (Maxwell RSC simplyRNA Cells Kit) från flera olika kända antal celler.

Cellberedning

Färskare celler skördades, tvättades med 1 × PBS och resuspenderades med 1 × 106 celler per ml i 1 × PBS och 0,04 % bovint serumalbumin. Färska PBMC:er frystes in vid 10 × genom att återuppta PBMC:er i DMEM+40% FBS+10% DMSO, frysas till -by °C i en CoolCell® FTS30 (BioCision) och sedan placeras i flytande kväve för förvaring.

Frusna cellflaskor från ALLCELLS och transplantationsstudier tinades snabbt upp i ett vattenbad vid 37 °C i ∼2 minuter. Flaskorna avlägsnades när en liten iskristall fanns kvar. Upptinade PBMC:er tvättades två gånger i mediet och resuspenderades sedan i 1 × PBS och 0,04 % bovint serumalbumin vid rumstemperatur. Cellerna centrifugerades vid 300 r.c.f. i 5 minuter varje gång. Upptinade BMMC:er tvättades och resuspenderades i 1 × PBS och 20 % FBS. Slutkoncentrationen av tinade celler var 1 × 106 celler per ml.

Konstruktion av sekvenseringsbibliotek med hjälp av GemCode-plattformen

Cellulära suspensioner laddades på ett GemCode Single-Cell Instrument (10x Genomics, Pleasanton, CA, USA) för att generera GEMs av enskilda celler. RNA-Seq-bibliotek för enskilda celler framställdes med GemCode Single-Cell 3′ Gel Bead and Library Kit (säljs nu som P/N 120230, 120231, 120232, 10x Genomics). GEM-RT utfördes i en C1000 Touch Thermal cycler med 96-Deep Well Reaction Module (Bio-Rad; P/N 1851197): 55 °C i 2 timmar, 85 °C i 5 minuter; hålls vid 4 °C. Efter RT bröts GEM:s och det enkelsträngade cDNA:t rengjordes med DynaBeads MyOne Silane Beads (Thermo Fisher Scientific; P/N 37002D) och SPRIselect Reagent Kit (0,6 × SPRI; Beckman Coulter; P/N B23318). cDNA:t amplifierades med hjälp av C1000 Touch Thermal cycler med 96-Deep Well Reaction Module: 98 °C i 3 minuter; 14 cykler: 98 °C i 15 s, 67 °C i 20 s och 72 °C i 1 minut; 72 °C i 1 minut; hålls vid 4 °C. Den amplifierade cDNA-produkten rengjordes med SPRIselect Reagent Kit (0,6 × SPRI). Därefter klipptes cDNA till ∼200 bp med hjälp av ett Covaris M220-system (Covaris; P/N 500295). Indexerade sekvenseringsbibliotek konstruerades med hjälp av reagenserna i GemCode Single-Cell 3′ Library Kit enligt följande steg: (1) reparation av ändar och A-tailing, (2) adapterligering, (3) rensning efter ligering med SPRIselect, (4) PCR och rensning av provindex. Stråkodssekvenseringsbiblioteken kvantifierades genom kvantitativ PCR (KAPA Biosystems Library Quantification Kit for Illumina platforms P/N KK4824). Sekvenseringsbiblioteken laddades vid 2,1 pM på en Illumina NextSeq500 med 2 × 75 paired-end-satser med följande läslängd: 98 bp Read1, 14 bp I7 Index, 8 bp I5 Index och 10 bp Read2. Vissa tidigare bibliotek gjordes med 5 nt UMI, och 5 bp Read2 erhölls i stället. Dessa bibliotek har dokumenterats i kompletterande tabell 1.

ERCC assay

ERCC syntetiska spike-in RNAs (Thermo Fisher Scientific; P/N 4456740) späddes ut (1:10 eller 1:50) och laddades i ett GemCode Single-Cell Instrument, som ersatte de celler som normalt används för att generera GEMs. Spike-in Mix1 och Mix2 testades båda. Ett något ändrat protokoll användes eftersom endast en liten del av GEMs samlades in för RT- och cDNA-amplifiering. Efter avslutad GEM-RT avlägsnades 1,25 μl av emulsionen och tillsattes till en bifasisk blandning av Recovery Agent (125 μl) (P/N 220016) och 25 mM additiv 1 (30 μl) (P/N 220074, 10x Genomics). Återvinningsmedlet avlägsnades sedan och den återstående vattenlösningen rengjordes med SPRISelect Reagent Kit (0,8 × SPRI). cDNA amplifierades med hjälp av C1000 Touch Thermal cycler med 96-Deep Well Reaction Module: 98 °C i 3 minuter; 14 cykler: 98 °C i 15 s, 67 °C i 20 s och 72 °C i 1 minut; 72 °C i 1 minut; hålls vid 4 °C. Amplifierad cDNA-produkt rengjordes med SPRIselect Reagent Kit (0,8 × ) cDNA klipptes därefter till ∼200 bp med hjälp av ett Covaris M220-system för att konstruera provindexerade bibliotek med 10x Genomics-adaptrar. Det förväntade antalet ERCC-molekyler beräknades utifrån den mängd ERCC-molekyler som användes och provutspädningsfaktorer. Antalet jämfördes med antalet detekterade molekyler (UMI-antal) för att beräkna konverteringseffektiviteten.

ddPCR-analys

Jurkat-celler användes i ddPCR-analyser för att uppskatta konverteringseffektiviteten enligt följande: (1) Mängden RNA per Jurkat-cell bestämdes genom att kvantifiera (Qubit, Invitrogen) RNA som extraherats (Maxwell RNA Purification Kits) från flera olika kända antal Jurkat-celler. (2) Bulk RT-ddPCR (Bio-Rad One-Step RT-ddPCR Advanced Kit for Probes 1864021) utfördes på det extraherade RNA:t för att bestämma antalet kopior per cell av åtta utvalda gener. (3) Ungefär 5 000 Jurkat-celler bearbetades med GemCode Single-Cell 3′-plattformen, och enkelsträngat cDNA samlades in efter RT i GEMs enligt de protokoll som anges i avsnittet ”Sequencing library construction using the GemCode platform”. cDNA-kopior av de åtta generna bestämdes med hjälp av ddPCR (Bio-Rad ddPCR Supermix for Probes (no dUTP) P/N 1863024). Det faktiska antalet Jurkat-celler fastställdes genom sekvensering av en delmängd av GEM-RT-reaktionerna på en MiSeq. Konverteringseffektiviteten är förhållandet mellan cDNA-kopior per cell (steg 3) och RNA-kopior per cell från RT-ddPCR i bulk (steg 2), med en effektivitet på 50 % i RT-ddPCR34.

Sondsekvenserna för ddPCR-analysen är följande: SERAC1_f, 5′-CACGAGAGCCGCCAGC-3′ och SERAC1_r, 5′-TCTGCAACACAGATGACGCAATAAG-3′; SERAC1_p: /56-FAM/CGCCCCTGCCG/ZEN/GCAGAATGTC/3IABkFQ/. AP1S3_f, 5′-GAAGCAGCCATGGTCTAAGC-3′ och AP1S3_r, 5′-CCTTGTCGACTGAAGAGCAATATATG-3′; AP1S3_p: /56-FAM/CGGCCCCCAGC/ZEN/CACGATGATGATACAT/3IABkFQ/OR. AOV1_f, 5′-CCGGAAGAGTGGGTCTCTCGTOR-3′ och AOV1_r, 5′-TTCTTCTCATAGCCTTCCCCCGATACCOR-3′; AOV1_p: /56-FAM/TCGGTGATGG/ZEN/CGGATGAGAGAGGTTTCA/3IABkFQ/. DOLPP1_f, 5′-ATGGCAGCGGACGGA-3′ och DOLPP1_r, 5′-GGCTCAGGTAGGCAAGGA-3′; DOLPP1_p: /56-FAM/CCACACGTCGA/ZEN/ATATATCCTGCAGGAGGTGATCT/3IABkFQ/. KPNA6_f, 5′-TGAAAGAGCTGCCGCTGAAG-3′ och KPNA6_r, 5′-CCCTGGGCTCGCCAT-3′; KPNA6_p: /56-FAM/CGGACCCCCGC/ZEN/GATGGAGACC/3IABkFQ/. ITSN2_f, 5′-GTGACACAGGCTACGCAACAG-3′ och ITSN2_r, 5′-TCCTGAGTTTTCCTCTTGCTAGCT-3′; ITSN2_p: /56-FAM/AGGGCGCCCA/ZEN/GATGGCTGA/3IABkFQ/. LCMT1_f, 5′-GTCGACCCCCCGCTTCCA-3′ och LCMT1_r, 5′-GGTCATGCCAGTAGCCAATG-3′; LCMT1_p: /56-FAM/ATGCTCTTCCCCC/ZEN/TGTGCAAGAGAGGTTTGTGC/3IABkFQ/. AP2M1_f, 5′-GCAGCGGCGGGCAGACG-3′ och AP2M1_r, 5′-ATGGCGGCAGCAGATCAGTCT-3′; AP2M1_p: /56-FAM/CATCGCCTCT/ZEN/GAGAACACAGACACCTGGTG/3IABkFQ/.

Beräkning av effektiviteten för cellfångst

Effektiviteten beräknas genom att ta förhållandet mellan antalet celler som detekteras genom sekvensering och antalet celler som laddas in i chipet. Det sistnämnda bestäms från (tillförd volym × inmatningskoncentration av celler). Ingångskoncentrationen av celler bestämdes med hjälp av en Countess II Automated Cell Counter (Thermo Fisher Scientific). Det är värt att notera att det finns ett fel på 15-20 % i cellräkningar, vilket kan förklara åtminstone en del av variationen i de beräknade verkningsgraderna.

Chimerism assay

PowerPlex 16 System (Promega) användes tillsammans med en Applied Biosystems (Life Technologies) 3130xl Genetic Analyzer. Donator BMMCs användes som referensbaslinje.

Allinering, streckkodstilldelning och UMI-räkning

Cell Ranger Single-Cell Software Suite användes för att utföra demultiplexering av prover, streckkodsbearbetning och 3′-genräkning av enskilda celler (http://software.10xgenomics.com/single-cell/overview/welcome). Först utfördes demultiplexering av provet baserat på den 8 bp avlästa provindexet för att generera FASTQ:er för Read1- och Read2-avläsningarna med parvis slutning samt GemCode-streckkoden på 14 bp. Tio baspars UMI-taggar extraherades från Read2 (14 bibliotek gjordes med 5 bp UMI-taggar, vilket noteras i den kompletterande tabellen 1, på grund av en tidigare iteration av metoderna. För dessa prover extraherades UMI-taggar på 5 bp från Read2). Därefter anpassades Read1, som innehåller cDNA-insatsen, till ett lämpligt referensgenom med hjälp av STAR35. För musceller användes mm10 och för mänskliga celler användes hg19. För prover med blandningar av mus- och människoceller användes föreningen av hg19 och mm10. För ERCC-prover användes ERCC-referensen (https://tools.thermofisher.com/content/sfs/manuals/cms_095047.txt).

Nästan filtrerades GemCode-streckkoder och UMIs. Alla kända förtecknade streckkoder som är 1-Hamming-avstånd från en observerad streckkod beaktas. Därefter beräknas den efterföljande sannolikheten för att den observerade streckkoden har producerats av ett sekvenseringsfel, med tanke på den observerade streckkodens baskvaliteter och den tidigare sannolikheten för att observera kandidatstreckkoden (hämtad från den totala fördelningen av antalet streckkoder). Om den efterföljande sannolikheten för någon kandidatstreckkod är minst 0,975 korrigeras streckkoden till den kandidatstreckkod som har den högsta efterföljande sannolikheten. Om alla kandidatsekvenser är lika sannolika, väljs den som visas först i lexikalisk ordning.

UMI:er med sekvenseringskvalitetspoäng >10 betraktades som giltiga om de inte var homopolymerer. Kval=10 innebär en noggrannhet på 90 % av basen. Ett UMI som är 1-Hamming-avstånd från ett annat UMI (med fler läsningar) för samma cellstreckkod och gen korrigeras till UMI med fler läsningar. Detta tillvägagångssätt är nästan identiskt med det i Jaitin et al.4 och liknar det i Klein et al.8 (även om Klein et al.8 också använde UMI:er för att lösa multimappade avläsningar, vilket inte genomfördes här).

Sist markerades PCR-dubbletter om två uppsättningar av avläsningspar delade en streckkodssekvens, en UMI-tagg och ett gen-ID (Ensembl GTF:erna GRCh37.82, ftp://ftp.ensembl.org/pub/grch37/release-84/gtf/homo_sapiens/Homo_sapiens.GRCh37.82.gtf.gz och GRCm38.84, ftp://ftp.ensembl.org/pub/release-84/gtf/mus_musculus/Mus_musculus.GRCm38.84.gtf.gz, användes). Endast säkert mappade (MAPQ=255), icke-PCR-duplikat med giltiga streckkoder och UMI användes för att generera gen-streckkodsmatris.

Streckkoder för celler bestämdes utifrån fördelningen av UMI-räkningar. Alla översta streckkoder inom samma storleksordning (>10 % av den översta n:e streckkoden, där n är 1 % av det förväntade återvunna cellantalet) betraktades som cellstreckkoder. Antal läsningar som ger meningsfull information beräknas som produkten av fyra mätvärden: (1) giltiga streckkoder, (2) giltigt UMI, (3) associerat med en cellstreckkod och (4) säkert mappat till exoner.

I blandningsexperimenten med mus och människa definierades multipletfrekvensen som två gånger frekvensen av cellstreckkoder med signifikanta UMI-räkningar från både mus och människa, där den översta 1 % av UMI-räkningarna ansågs vara signifikanta. Omfattningen av barcode crosstalk bedömdes av fraktionen musläsningar i mänskliga barcodes, eller vice versa.

Prover som bearbetas från flera kanaler kan kombineras genom att sammanfoga gen-cell-barcode-matriser. Denna funktionalitet finns i Cell Ranger R Kit (http://support.10xgenomics.com/single-cell/software/pipelines/latest/rkit). Sekvenseringsdata från flera sekvenseringskörningar av ett bibliotek kan kombineras genom att räkna icke-duplicerade läsningar. Denna funktion finns i Cell Ranger pipeline. Dessutom kan sekvenseringsdata underprovtas för att erhålla ett givet antal UMI-räkningar per cell. Denna funktionalitet finns också i Cell Ranger R Kit och är användbar när man kombinerar data från flera prover för jämförelse.

PCA-analys av blandning av Jurkat- och 293T-celler

Gen-cell-streckkodsmatris från vart och ett av de fyra proverna sammanfogades. Endast gener med minst en UMI-räkning som upptäckts i minst en cell används. UMI-normalisering utfördes genom att först dividera UMI-räkningarna med de totala UMI-räkningarna i varje cell, följt av en multiplikation med medianen av de totala UMI-räkningarna över cellerna. Därefter tog vi den naturliga loggen av UMI-räkningarna. Slutligen normaliserades varje gen så att medelsignalen för varje gen är 0 och standardavvikelsen 1. PCA kördes på den normaliserade matrisen av gener och streckkoder. De normaliserade UMI-räkningarna för varje gen används för att visa uttrycket av en markör i en tSNE-plott.

SNV-analys av Jurkat- och 293T scRNA-seq-data

SNV:er kallades genom att köra Freebayes 1.0.2 (ref. 36) på den genome BAM som producerats av Cell Ranger. SNV av hög kvalitet (SNV calling Qual>=100 med minst 10 UMI-räkningar från minst två celler; indels ignoreras) som endast observerades i Jurkat- eller 293T-celler (men inte i båda) valdes ut. Cellerna märktes som Jurkat eller 293T baserat på Jurkat- och 293T-specifika SNV-räkningar, där andelen räkningar från den andra arten är <0,2. Celler med en fraktion av SNV från någon av arterna mellan 0,2 och 0,8 betraktas som multipletter. Den härledda multipletfrekvensen är 2* den observerade multipletfrekvensen (för att ta hänsyn till Jurkat:Jurkat och 293T:293T multiplets).

PCA- och tSNE-analys av PBMCs

Gener med minst ett UMI-räkningstal som detekterats i minst en cell används. De 1 000 mest variabla generna identifierades utifrån deras medelvärde och spridning (varians/medelvärde), vilket liknar det tillvägagångssätt som används av Macoscko et al.7 Generna placerades i 20 bins utifrån deras medeluttryck. Den normaliserade spridningen beräknas som den absoluta skillnaden mellan spridningen och medianens spridning av uttrycksgenomsnittet, normaliserad med medianens absoluta avvikelse inom varje bin.

PCA kördes på den normaliserade gen-barkods-matrisen för de 1 000 mest variabla generna för att minska antalet dimensioner för funktionen (genen). UMI-normalisering utfördes genom att först dividera UMI-räkningar med de totala UMI-räkningarna i varje cell, följt av multiplikation med medianen av de totala UMI-räkningarna över cellerna. Därefter tog vi den naturliga logaritmen av UMI-räkningarna. Slutligen normaliserades varje gen så att medelsignalen för varje gen är 0 och standardavvikelsen 1. PCA kördes på den normaliserade matrisen av gener och streckkoder. Efter att ha kört PCA utfördes Barnes-hut37 approximation till t-SNE16 på de första 50 PC:erna för att visualisera cellerna i ett tvådimensionellt utrymme. Femtio PCs användes eftersom: (1) att använda alla PCs skulle ta mycket lång tid med tSNE-analys, (2) de förklarade ∼25 % av den totala variansen. K-means15-klustering kördes för att gruppera celler för klusteranalysen. k=10 valdes baserat på summan av kvadratfelet scree plot (Supplementary Fig. 5d).

Identifiering av klusterspecifika gener och markörbaserad klassificering

För att identifiera gener som är berikade i ett specifikt kluster beräknades medeluttrycket av varje gen över alla celler i klustret. Därefter jämfördes varje gen från klustret med medianuttrycket av samma gen från celler i alla andra kluster. Generna rangordnades utifrån skillnaden i uttryck och de tio mest berikade generna från varje kluster valdes ut. För hierarkisk klustring beräknades parvis korrelation mellan varje kluster, och centrerat uttryck för varje gen användes för visualisering med hjälp av värmekarta.

Klassificering av PBMCs härleddes från annoteringen av klusterspecifika gener. I fallet med kluster 10 upptäcktes marköruttryck av flera celltyper (t.ex. B, dendritiska och T). Eftersom den relativa klusterstorleken för B, dendritisk och T är 5,7 %, 6,6 % respektive 81 %, skulle vi förvänta oss att kluster 10 (som bara är 0,5 %) skulle innehålla multipletter bestående främst av B:dendritisk (0.36 %) och B:dendritic:T (0,3 %).

Selektion av renade subpopulationer av PBMCs

Varje population av renade PBMCs downsamplades till ∼16k reads per cell. PCA, tSNE och k-means-klustering utfördes för varje nedskalad matris, enligt samma steg som beskrivs i PCA- och t-SNE-analys av PBMCs. Endast ett kluster upptäcktes i de flesta prover, vilket överensstämmer med FACS-analyserna (kompletterande figur 6). För prover med mer än ett kluster valdes endast de kluster som uppvisade det förväntade markörgenuttrycket ut för nedströmsanalys. För CD14+ monocyter observerades två kluster som identifierades som CD14+ monocyter och dendritiska celler baserat på uttryck av markörgenerna FTL respektive CLEC9A.

Cellklassificeringsanalys med hjälp av renade PBMC:er

Varje population av renade PBMC:er downsamplades till ∼16k säkert mappade reads per cell. Därefter beräknades en genomsnittlig (medel) genuttrycksprofil över alla celler. Därefter jämfördes genuttrycket från varje cell i den komplexa populationen med genuttrycksprofilerna för renade populationer av PBMCs genom Spearmans korrelation. Cellen tilldelades ID för den renade populationen om den hade den högsta korrelationen med den populationen. Observera att skillnaden mellan den högsta och näst högsta korrelationen var liten för vissa celler (t.ex. skillnaden mellan cytotoxiska T- och NK-celler), vilket tyder på att celltilldelningen inte var lika säker för dessa celler. Några av de renade PBMC-populationerna överlappade varandra. Till exempel omfattar CD4+ T-hjälparceller alla CD4+ celler. Detta innebär att celler från detta prov kommer att överlappa med celler från prov som innehåller CD4+-celler, inklusive CD4+/CD25+ T reg, CD4+/CD45RO+ T memory, CD4+/CD45RA+/CD25- naive T. När en cell tilldelades ID:t CD4+ T-hjälparcell på grundval av korrelationspoängen, kontrollerades alltså den näst högsta korrelationen för att se om det var ett av CD4+-proven. Om så var fallet uppdaterades cellens ID till celltypen med nästa högsta korrelation. Samma procedur utfördes för CD8+ cytotoxisk T och CD8+/CD45RA+ naiv cytotoxisk T (som är en undergrupp av CD8+ cytotoxisk T).

Den R-kod som användes för att analysera 68k PBMCs och renade PBMCs finns här: https://github.com/10XGenomics/single-cell-3prime-paper.

Cellklustring och klassificering med Seurat

Gen-cell-barcode-matrisen för 68k PBMCs log-transformerades som indata till Seurat. De 469 mest variabla generna som valdes ut av Seurat användes för att beräkna PCs. De första 22 PCs var signifikanta (P<0,01) baserat på den inbyggda jackstraw-analysen och användes för tSNE-visualisering. Cellklassificering togs från Cellklassificeringsanalys med hjälp av renade PBMCs.

Variation mellan färska och frysta PBMCs

Sekvenseringsdata från 68k färska PBMCs och 3k frysta PBMCs nedprovades så att varje prov har ∼14k säkert mappade läsningar per cell. Endast gener som detekteras i minst en cell inkluderades i jämförelsen, som använder medelvärdet för varje gen över alla celler.

För jämförelsen av cellklassificering mellan renade och frysta PBMC:er, sammanförde vi alla celler som var märkta som T- eller natural killer-celler tillsammans. Detta beror på att subpopulationerna inom T och mellan T och naturliga mördarceller ibland är svåra att klustra separat. Vi ville inte att jämförelsen mellan färska och frysta celler skulle påverkas av de klustermetoder som användes.

SNV-baserad genotyptilldelning

SNV:er kallades genom att köra Freebayes 1.0.2 (ref. 36) på den genom BAM som producerades av Cell Ranger. Endast SNV:er med stöd från minst två cellstreckkoder, med en minimal SNV Qual-poäng >=30, minimal SNV base Qual>=1 inkluderades. Antalet referensalleler (R) och alternativa alleler (A) beräknades för varje SNV, vilket gav en matris med cellreferens-UMI-räkningar och cell-alternativ-allel-UMI-räkningar. Dessa matriser modellerades som en blandning av två genomer där sannolikheten för någon av de tre genotyperna (R/R, R/A eller A/A) på en plats antogs vara binomiskt fördelad med en fast felprocent på 0,1 %. För varje prov har två modeller tagits fram parallellt, en där endast ett genom finns (K=1) och en där två genom finns (K=2). Inferens av modellparametrarna (cell-till-genom-tilldelning och K uppsättningar genotyper) utfördes med hjälp av en Gibbs sampler för att approximera deras posteriora fördelningar. För att förbättra problemet med etikettbyte vid Monte Carlo-inferens av blandningsmodeller utfördes en ny etikettering av de samplade cell-till-genom-tilldelningarna i enlighet med Stephens et al.38

I in silico-cellblandningsexperiment, när K=2-modellen misslyckades med att på ett adekvat sätt separera de två genomerna, rapporterade den en fördelning av sannolikheter i efterhand som låg nära 0,5 för cell-till-genom-tilldelningarna, vilket tyder på bristande förtroende för dessa tilldelningar. Vi tillämpade ett krav på att 90 % av cellerna ska ha en efterföljande sannolikhet >75 % för att välja K=2-modellen framför K=1-modellen. Att välja K=1 indikerar att blandningsfraktionen ligger under metodens detektionsnivå, som i in silico-blandningsexperiment fastställdes till 4 % av 6 000 celler.

Genotypjämförelse med det rena provet

För att säkerställa tilldelningen av genotyper till individer beaktades endast delade SNV:er mellan genotypgruppen och det rena provet. Därefter jämfördes den genomsnittliga genotypen för alla celler med genotypen för det rena provet. För att få en viss baslinje för den procentuella genotypöverlappningen mellan olika individer utförde vi parvisa jämförelser av genotyper som kallats från samma individer (11 parvisa jämförelser) eller från olika individer (15 parvisa jämförelser). Den procentuella genotypöverlappningen mellan samma individer är i genomsnitt ∼98±0,3 %, medan den procentuella genotypöverlappningen mellan olika individer är i genomsnitt ∼73±2 %.

PCA- och tSNE-analys av BMMCs

Data från sex prover användes: två friska kontroller, AML027 före och efter transplantation samt AML035 före och efter transplantation. Varje prov nedprovades till ∼10k säkert mappade läsningar per cell. Därefter sammanfogades matrisen med streckkoder för gener och celler från varje prov. PCA, tSNE och k-means-klustering utfördes på den sammanslagna matrisen, enligt samma steg som beskrivs i PCA- och tSNE-analysen av PBMCs. För k-means-klustering användes K=10 baserat på krökningen i summan av kvadratfelet scree plot.

Klusterspecifika gener identifierades enligt de steg som beskrivs i ”Identifiering av klusterspecifika gener och markörbaserad klassificering”. Klassificering tilldelades baserat på klusterspecifika gener och baserat på uttryck av några välkända markörer för immuncellstyper. ’Blasts and Immature Ery 1’ avser kluster 4, som uttrycker CD34, en markör för hematopoietiska progenitorer39 , och Gata2, en markör för tidiga erytroider40. ”Immature Ery 2” avser kluster 5 och 8, som visar uttryck av Gata1, en transkriptionsfaktor som är nödvändig för erytropoesin41 , men inte CD71, som ofta återfinns i mer etablerade erytroida celler39. ”Immature Ery 3” avser kluster 1, som visar uttryck av CD71. ”Mature Ery” avser kluster 2. HBA1, en markör för mogna erytroida celler, upptäcks företrädesvis i kluster 2. Kluster 3 tilldelades beteckningen ”omogna granulocyter” på grund av uttrycket av tidiga granulocytmarkörer som AZU1 och IL8 (ref. 42) och avsaknaden av uttryck av CD16. Kluster 7 klassificerades som ”monocyter” på grund av uttrycket av t.ex. CD14 och FCN1. ”B” avser kluster 6 och 9 på grund av markörer som CD19 och CD79A. ’T’ hänvisar till kluster 10 på grund av markörer som CD3D och CD8A.

Datatillgänglighet

Alla relevanta data är tillgängliga från författarna. RNA-seq-data från enskilda celler har deponerats i Short Read Archive under accession nummer SRP073767. Data finns också tillgängliga på http://support.10xgenomics.com/single-cell/datasets. Analyskoden för 68k PBMC-analysen finns tillgänglig på https://github.com/10XGenomics/single-cell-3prime-paper.

.