Geni umani che sfuggono all’inattivazione X rivelati dai dati di espressione delle singole cellule

Un quadro per misurare la fuga dall’inattivazione X in singole cellule

Abbiamo identificato gli sfuggenti analizzando l’espressione genica da singole cellule somatiche usando la metodologia scRNA-Seq (vedi metodi). Per valutare la sensibilità del metodo, confrontiamo l’espressione del cromosoma X (ChrX) con altri cromosomi autosomici. In particolare, ci siamo concentrati sul cromosoma ricco di geni 17 (Chr17) come prototipo di un cromosoma autosomico. Chr17 è stato selezionato in quanto rappresenta un cromosoma con un numero minimo di geni impressi specifici dei genitori. Le proprietà quantitative di ChrX e Chr17 sono elencate in Fig. 1a.

Questo studio si basa sull’analisi di due risorse di origine femminile: (i) Fibroblasti primari UCF1014 (con 104 cellule, vedi Metodi). Questo set è specificato da una maggiore copertura dei dati trascrittomici, ma manca di informazioni sulla fasatura degli aplotipi (Fig. 1b); (ii) un set di dati più piccolo di linfoblasti clonali (n = 25) dalla linea cellulare GM12878 con genomi diploidi parentali completamente fasati e sequenziati (Fig. 1c). In entrambi i set di dati, la trascrizione a SNPs eterozigoti (hSNPs) è la fonte di informazioni per determinare l’espressione monoallelica o biallelica. Ogni hSNP, in ogni cella, che è supportato da prove di espressione sopra una soglia predeterminata è considerato un SNP informativo (iSNP) (vedi Metodi, file aggiuntivo 1: testo). La somma di iSNPs per gene definisce la sua etichetta unica come gene inattivato o escapee (vedi Metodi, Fig. 1b-c, file aggiuntivo 1: testo).

Quantificare l’espressione biallelica da fibroblasti primari a singola cellula

Abbiamo analizzato i dati scRNA-Seq pubblicati da fibroblasti umani primari femminili. Prima dell’analisi, ci siamo occupati di un trabocchetto sperimentale rilevante per molti studi su singole cellule. Il trabocchetto riguarda le doppie cellule in cui più di una cella è sequenziata per una libreria. In un tale scenario, diversi cromosomi X attivi (Xa) da due diverse cellule che sono inclusi nel campione produrrà un segnale biallelico lungo l’intero cromosoma X. Anche se la frazione di doppietti dovrebbe essere piccola, può portare a un’interpretazione errata. Pertanto, prima di analizzare i dati dai fibroblasti, abbiamo rivisitato tutti i 104 fibroblasti e testato il loro rapporto biallelico rispetto a ChrX (vedi Metodi). Tre cellule hanno mostrato un grado eccezionalmente alto di espressione biallelica che potrebbe indicare una miscela di due cromosomi X parentali (file aggiuntivo 1: testo e figura S3). Abbiamo rimosso tutte e tre le cellule sospette da tutte le analisi.

Poi, per ogni singola cella, abbiamo contato il numero di letture che sono state univocamente mappate agli alleli hSNP. Il rapporto allelico (AR) per ogni iSNP è definito come la frazione delle letture mappate all’allele alternativo (Alt) sul totale delle letture (vedi metodi, file aggiuntivo 2: tabella S1). La figura 2a-c riassume l’AR di ChrX, Chr17, e gli interi cromosomi autosomici secondo la collezione di fibroblasti primari (101 su 104 cellule). Inoltre, la Fig. 2d mostra la distribuzione dell’AR di un set annotato di geni imprinted dai tessuti della pelle (secondo ). Come precedentemente riportato, una distorsione nella mappatura verso il genoma di riferimento (AR = 0) è evidente (Fig. 2a-d). Inoltre, una frazione sostanziale di espressione monoallelica è stata osservata per tutti i set testati (Fig. 2a-d). Questa apparizione dominante di espressione monoallelica in cellule singole è causata da una combinazione sia di sotto-campionamento di trascrizioni, sia di un fenomeno noto come “bursting trascrizionale”.

Ci siamo concentrati solo su iSNPs che mostrano una firma non monoallelica (cioè, escludendo AR = 0 e AR = 1). Abbiamo osservato una marcata differenza nella distribuzione AR di ChrX e dei geni imprinted rispetto a Chr17 e tutti i cromosomi autosomici (confrontare Fig. 2a a b e Fig. 2c a d). Di conseguenza, diverse osservazioni dai risultati mostrati in Fig. 2 possono essere tratte: (i) Chr17 e tutti gli autosomi condividono un profilo AR simile. (ii) Una chiara tendenza verso l’espressione equilibrata (AR = 0.5) è evidente per tutti i cromosomi autosomici (Fig. 2a-b), ma non ChrX o geni imprinted (Fig. 2c-d). (iii) La frazione di espressione non-monoallelica nei cromosomi autosomici è più alta (~ 18%) rispetto a ChrX (~ 9%). (iv) La frazione di espressione non monoallelica nei geni imprinted mostra un livello intermedio (13%). Tale livello intermedio è probabilmente un riflesso dell’incoerenza inerente all’identità dei geni imprinted. File aggiuntivo 3: La tabella S2 elenca gli iSNPs di supporto per tutti i cromosomi analizzati nei fibroblasti compreso il set di geni imprinted.

Identificare gli escapees nei fibroblasti primari a cellula singola

Nel dataset dei fibroblasti primari, ci sono 232 e 485 geni che sono supportati da prove iSNPs per ChrX e Chr17, rispettivamente. Poiché queste cellule mancano di informazioni sulla fasatura del genoma (Fig. 1b), le informazioni sulla fuga dal cromosoma X sono limitate al set di iSNPs biallelici (vedi Fig. 1b). Abbiamo aggregato gli iSNPs secondo i loro geni corrispondenti (Fig. 1b). L’aggregazione viene eseguita su diverse cellule singole e su più iSNPs all’interno di una specifica coppia cellula-gene. Un gene sarà etichettato candidato escapee quando è associato a più iSNPs biallelici. Complessivamente abbiamo identificato 24 geni di questo tipo (Tabella 1) che rappresentano il 10,3% di tutti i geni espressi in ChrX. Come previsto, la frazione di geni su Chr17 che mostra l’espressione biallelica è sostanzialmente più alta (49,3%, file aggiuntivo 4: Tabella S3).

Tabella 1 elenca i candidati alla fuga a vari gradi di supporto. Per esempio, i geni ZFX (Zinc finger X-chromosomal protein) e SMC1A (Structural maintenance of chromosomes protein 1A) sono fortemente supportati con 103 e 19 iSNPs biallelici, rispettivamente. Un ulteriore aumento dell’affidabilità dell’identificazione degli escape si basa sull’avere almeno 2 cellule indipendenti che hanno contribuito con informazioni sull’espressione biallelica. Mostriamo che 21 dei 24 geni hanno soddisfatto questo criterio rigoroso (Tabella 1). In particolare, tra gli escapees identificati abbiamo rilevato solo tre geni PAR (SLC25A6, CD99 e DHRSX, Tabella 1). L’assegnazione di questi geni come escape concorda con l’espressione PAR prevista. Dal numero di geni PAR biallelici sui geni PAR espressi, abbiamo stimato il tasso di scoperta falso negativo per gli escapees per essere alto come 70% (cioè, mancato 7 di 10 geni PAR espressi). File aggiuntivo 4: La tabella S3 mostra il supporto per la tabella 1.

Quantificare l’espressione allelica da linfoblasti clonali fasati

Un limite importante nel protocollo descritto sopra riguarda la mancanza di fasatura dell’aplotipo dei genitori. In questa impostazione, iSNPs non può essere assegnato a Xa o Xi. Di conseguenza, la scelta casuale di Xi che caratterizza le cellule primarie limita il tasso di scoperta per gli escapees. Abbiamo esteso l’analisi di scRNA-Seq ai linfoblasti di origine femminile della linea cellulare clonale GM12878 . Vedi Fig. 1c e file aggiuntivo 1: Figura S4.

La figura 3a mostra il profilo di espressione per cellule singole di linfoblasti clonali (n = 25) (file aggiuntivo 2: Tabella S1, file aggiuntivo 5: Tabella S4). In ogni singola cellula, l’espressione monoallelica riflette la combinazione di un sottocampionamento di trascrizioni e il fenomeno del bursting trascrizionale. È chiaro che l’espressione materna dalla Xa domina (Fig. 3a, in alto). Un’osservazione che concorda con l’origine materna Xa riportata della linea cellulare GM12878 . Nella maggior parte delle cellule analizzate, una piccola ma sostanziale frazione dell’espressione totale misurata proviene dal cromosoma Xi paterno (Fig. 3a, in alto). Al contrario, Chr17 e i cromosomi autosomici mostrano un’espressione uguale da entrambi gli alleli (Fig. 3a, pannelli centrale e inferiore).

Dalla figura, è evidente che il fenomeno del burst trascrizionale interessa tutti i cromosomi compreso ChrX. Per valutare l’impatto di questo fenomeno sull’identificazione dei geni come escapees, abbiamo confrontato singole cellule rispetto a un pool di cellule (Pool100, Fig. 3a, barra di destra). Mentre la maggior parte degli iSNPs di Chr17 mostrano profili biallelici, ChrX rimane dominata da un’espressione monoallelica materna.

La figura 3b è una visione aggregata di ChrX, Chr17 e dei cromosomi autosomici. I dati si basano su 375, 808 e 20.212 hSNPs espressi, rispettivamente. La figura 3b (pannelli centrale e destro) mostra una partizione uguale degli alleli parentali di Chr17 e dei cromosomi autosomici (Fig. 3b, in alto). Eseguendo la stessa analisi sui dati raccolti da Pool100 (Fig. 3c) mostra che la partizione degli alleli parentali rimane praticamente invariata (confrontare la frazione occupata dai colori rosa e blu, Fig. 3b-c). Inoltre, abbiamo osservato uno spostamento da un’espressione monoallelica (Fig. 3b-c, colore pieno) a una biallelica (Fig. 3b-c, colore spogliato). La frazione dell’espressione biallelica per Ch17 è aumentata dal 19% in singole cellule all’80% in Pool100, e per i cromosomi autosomici dal 18 al 79% (Fig. 3c, pannelli centrale e destro). I risultati da Pool100 indicano che l’espressione monoallelica osservata in singole cellule è praticamente abolita dalla media del segnale.

I risultati da ChrX (Fig. 3b (sinistra) sono fondamentalmente diversi rispetto a Chr17 o i cromosomi autosomici (Fig. 3b-c). La differenza più notevole è che solo il 21% degli iSNPs espressi sono associati all’allele Xi paterno in ChrX (Fig. 3b, in alto). Inoltre, nell’analizzare Pool100, la frazione di espressione biallelica rimane limitata (un passaggio dal 9% nelle singole cellule al 34% in Pool100). Il modello osservato di ChrX da Pool100 (Fig. 3c, a sinistra) è meglio spiegato da una media del segnale stocastico monoallelico (allo stesso grado degli altri cromosomi) pur mantenendo un forte segnale dell’espressione monoallelica Xa. Si veda il file aggiuntivo 5: Tabella S4 per il rapporto allelico dei linfoblasti di tutti i cromosomi testati e il Pool100.

Identificare le fughe dai linfoblasti monocellulari

La figura 4a è una vista gene-centrica che mostra la partizione allelica iSNP dei linfoblasti (colorata secondo la loro origine come espressione materna, paterna o mista, vedi Metodi). Solo il sottoinsieme dei geni che sono supportati da più iSNPs è elencato secondo il loro ordine lungo i cromosomi. Complessivamente riportiamo 93 geni annotati su ChrX (Fig. 4a, 30 escapees e 63 geni inattivati). Si noti che i geni X-inattivati rappresentano i geni che sono espressi principalmente dalla Xa materna. Un cluster di geni con un’espressione paterna all’estremità del braccio p di ChrX rappresenta l’espressione biallelica attesa dai geni PAR (Fig. 4a). Ulteriori prove di espressione paterna sono localizzate al XIC con geni come XIST, JPX e FTX. Mentre la maggior parte degli escape sono supportati da un numero limitato di iSNPs, alcuni di loro come ZFX, CD99, e SLC25A6 sono supportati da un numero relativamente elevato di iSNPs di supporto (48, 38 e 34, rispettivamente).

Un metodo alternativo per valutare l’entità del fenomeno di X-inattivazione è quantificare le prove direttamente dalla somma di tutte le letture sequenziate (abbreviato come il protocollo basato sulla lettura). La figura 4b-c confronta i conteggi delle letture da Chr17 (Fig. 4b) e ChrX (Fig. 4c) in base all’origine paterna rispetto a quella materna. Confrontiamo i dati di espressione dalle singole cellule e dal Pool100. La regressione lineare per l’espressione genica da Chr17 mostra un’alta correlazione fit-line (r2 = 0,823, Fig. 4b). Come previsto, la correlazione è più forte nei dati originati dal Pool100 (r2 = 0,946, Fig. 4b). Abbiamo concluso che nonostante l’espressione monoallelica dovuta al fenomeno del bursting trascrizionale, l’espressione allelica bilanciata di tutti i geni è fortemente sostenuta. Per ChrX tuttavia, la regressione lineare risultante delle singole cellule è scarsa (r2 = 0,238, Fig. 4c), e non è stata migliorata dai dati del Pool100 (r2 = 0,222, Fig. 4d). Ispezionando i dati di espressione per ChrX si nota che le linee di regressione si inclinano effettivamente verso l’espressione Xa materna (asse delle x). I dati di espressione sono coerenti con due linee di regressione distinte per ChrX. Uno che corrisponde ai geni inattivati (parallelo all’asse x), e l’altro corrisponde un’espressione biallelica affidabile.

Applicando il protocollo conservativo basato su iSNP porta all’identificazione di 30 geni come candidati di fuga che sono anche supportati dal protocollo basato sulla lettura (Tabella 2). Il read-based (cioè, etichettando un gene come fuggitivo per avere un numero minimo di letture paterne, vedi Metodi) ha ampliato la lista dei candidati fuggitivi per includere complessivamente 49 geni (file aggiuntivo 4: Tabella S3).

Testare l’origine parentale degli alleli lungo un gene nella stessa cella è un test rigoroso per l’affidabilità degli iSNPs. Questo test è valido solo per i geni con più hSNPs. Tali geni che sono supportati con due o più hSNPs esprimenti rappresentano il 44% dei geni. Consideriamo un gene coerente se l’espressione lungo il gene in una cella specifica non è monoallelica per entrambi gli alleli. Complessivamente, abbiamo identificato 3 geni incoerenti – TEX11, FTX e ZCCHC16. Per altri 6 geni, l’incoerenza è solo parziale poiché ci sono altre osservazioni di espressione biallelica. La stima dell’incoerenza completa (3 dei 29 geni che erano ammissibili per questo test) suggerisce un limite superiore per un’interpretazione errata del 10%. Inoltre, il protocollo basato su iSNP ha identificato 9 degli 11 geni PAR espressi. Quindi, estrapoliamo il tasso di rilevamento degli evasi all’82%. È interessante notare che l’analisi di Chr17, sotto l’ipotesi che non ci sia un bias allelico sistematico, ha mostrato che il 7,3 e il 9,6% dei geni erano associati all’espressione monoallelica materna e paterna, rispettivamente. Questi risultati forniscono un limite superiore del 17,9% alla probabilità di falsa etichettatura dei geni in Chr17 e possono essere utilizzati per stimare la limitazione del metodo.

Richiedere l’evidenza paterna da almeno due cellule riduce il numero di fughe da 49 a 18 (compresi 5 geni PAR). Molti geni a prova singola derivanti da una cella outlier specifica (SRR764803) che ha contribuito a 17 dei 49 geni riportati (Tabella 2, file aggiuntivo 5: Tabella S4). Questa cellula mostra i più alti livelli di espressione sia su ChrX che su Chr17 (file aggiuntivo 1: testo e figura S4). La tabella 2 elenca i candidati di fuga scoperti insieme alle loro prove di supporto (Tabella 2).

Confronto delle fughe identificate con le conoscenze attuali

Abbiamo testato la corrispondenza tra le fughe identificate dal nostro studio e un catalogo basato sulla letteratura. Questo catalogo unificato è stato compilato dall’integrazione di quattro studi indipendenti che copre 1144 geni da ChrX. I geni in questo catalogo sono suddivisi manualmente in nove categorie definite (vedi Metodi). La più grande rappresenta i geni che mancano di informazioni (45%). Circa il 15% dei geni (168/1144) sono considerati ‘associati alla fuga’ (vedi Metodi). Consideriamo l’insieme compilato di geni ‘associati agli escapee’ come un gold standard per testare il tasso di scoperta degli escapee nel nostro studio (totale di 124 geni, esclusi i geni PAR, chiamati collettivamente Balaton-Esc).

Abbiamo applicato un test statistico ipergeometrico (vedi Metodi) per valutare la sovrapposizione delle diverse liste di geni escapee (Fig. 5). La figura 5a mostra il numero di fughe identificate nei fibroblasti e nei linfoblasti (esclusi i geni PAR). Si noti che solo i geni che sono inclusi nel benchmark Balaton-Esc sono inclusi in questa analisi (Fig. 5). La Figura 5b mostra la significatività statistica della sovrapposizione tra le liste di geni della Fig. 5a e il Balaton-Esc . Come si può vedere, c’è una sovrapposizione significativa tra le fughe dai linfoblasti (Tabella 2 e File aggiuntivo 4: Tabella S3) e la lista Balaton Esc (Fig. 5b, p-value = 7.43E-8). Applicando lo stesso test per i fibroblasti primari (Tabella 1) ha portato ad una minore significatività (p-value = 4.07E-2).

La figura 5c mostra i geni di sovrapposizione tra gli escape identificati nel nostro studio e il catalogo Balaton-Esc (168 geni compresi i geni PAR). Abbiamo anche incluso una risorsa complementare basata su 940 trascrittomi da scRNA-Seq (25 geni di fuga, Tukiainen-Esc) . Il diagramma di Venn mostra che ognuno degli studi di cui sopra contribuisce alla conoscenza attuale sulle fughe. Le fughe dalle due risorse esterne si sovrappongono per 18 dei 25 geni riportati (72%). Come mostrato in Fig. 5c, il 62% delle fughe riportate dai linfoblasti si sovrappongono alle liste di fughe esterne, mentre i fibroblasti sono supportati solo dal 38% di sovrapposizione. In particolare, la maggior parte dei nostri geni di fuga candidati scoperti dai fibroblasti (62%) non hanno alcuna corrispondenza con le altre liste testate (file aggiuntivo 6: Tabella S5).

LncRNA estendono la lista dei candidati di fuga

Abbiamo esteso la nostra analisi agli RNA lunghi non codificanti (lncRNA) utilizzando gli stessi criteri di fuga utilizzati dai geni codificanti (Fig. 6). Complessivamente abbiamo identificato 15 lncRNA come candidati di fuga, tra i quali solo alcuni sono stati precedentemente studiati. La posizione dei lncRNA e delle fughe codificanti lungo ChrX è mostrata (Fig. 6a). Abbiamo testato le posizioni delle fughe lungo la ChrX rispetto a tutti i geni ChrX. Mentre la distribuzione posizionale per lncRNA fughe è simile (Kolmogorov-Smirnov test, p-value = 0.57), è diverso per codifica fughe (Kolmogorov-Smirnov test, p-value = 0.004, Fig. 6a).

Figura 6b elenca tutti i 15 geni lncRNA identificati escapee, tra cui, sono geni ncRNA dal XIC che coordinano l’attivazione e il mantenimento di X-inattivazione. Molti degli lncRNA sono localizzati in segmenti trascrizionalmente attivi (ad esempio, all’interno della PAR o della XIC), mentre altri sono localizzati in regioni non conservate che sono arricchite con ncRNA lunghi e corti. Per ulteriori lncRNA, compresi i geni inattivati si veda il file aggiuntivo 7: Tabella S6.

Partizione basata sull’evidenza dei geni fuggitivi

La figura 6c riassume la partizione secondo l’evidenza per tutti i geni che sono riportati con qualsiasi livello di fiducia come fuggitivi (file aggiuntivo 6: Tabella S5). Questo elenco comprende 75 candidati che sono riportati in questo studio, tra cui prove da Pool100, una raccolta di lncRNAs romanzo. Figura 6d fornisce un elenco dettagliato della scoperta da questi dati in vista di serval risorse esterne. Il numero di cellule che forniscono prove è anche indicato per ogni gene.

Tenendo conto di queste risorse esterne, abbiamo abbinato ogni gene secondo la qualità del supporto indipendente associato ad esso (Fig. 6c e d). In particolare, abbiamo suddiviso la lista dei 75 geni in quattro gruppi: (i) I geni sono etichettati come “confermati” se sono riportati come fughe da entrambe le risorse esterne precedentemente discusse . Ci sono 31 geni di questo tipo, per i quali questo studio fornisce un’ulteriore conferma della loro identità come escapees. (ii) Altri 7 geni sono etichettati come “approvati”. Questi geni sono etichettati come escapees da una sola delle due risorse esterne. Per questi geni, le prove indipendenti di questo studio approvano la loro identità. (iii) Altri 12 geni sono contrassegnati come fughe “probabili”. Questi geni sono assegnati in base all’accordo con un ulteriore rapporto esterno che riporta 114 fughe, che non è stato incluso nel catalogo basato sulla letteratura, e quindi può essere considerato come una risorsa indipendente. (iv) Altri 25 geni sono contrassegnati come fughe “possibili”. Questi geni mancano di qualsiasi prova di letteratura per sostenere la loro identità, quindi la loro assegnazione come fuggitivi rimane meno supportata. Tra questi, 5 geni sono supportati solo dalla cellula outlier del linfoblasto e sono quindi probabilmente falsi. In questo set, ci sono altri 10 geni che non sono stati riportati da nessuna delle tre risorse esterne discusse, la maggior parte dei lncRNA trascurati appartengono a questo gruppo. Due di questi lncRNA trascurati (TCONS_00017125, e TCONS_00017281, Fig. 6b) si trovano nella regione PAR che è eccezionalmente attiva nella trascrizione biallelica. I dati di origine per la Fig. 6d sono nel file aggiuntivo 6: Tabella S5.

.