ChIP-chip versus ChIP-seq: Lezioni per il disegno sperimentale e l’analisi dei dati

Analisi dei profili del DNA in ingresso

Per capire le differenze tecnologiche tra i dati ChIP basati sul sequenziamento e quelli basati sul microarray, abbiamo prima analizzato i profili dei frammenti di DNA reticolati e sonicati (DNA in ingresso) generati dal microarray (INPUT-chip) e dal sequenziamento ad alto rendimento (INPUT-seq). Poiché il profilo del DNA di input dovrebbe essere indipendente dall’anticorpo utilizzato per ChIP, questo confronto può dare un’idea delle differenze specifiche tra queste due tecnologie di profilazione. Abbiamo ottenuto INPUT-chip dati dal canale di fondo dei nostri dati microarray a due canali. Mentre questa piattaforma microarray utilizza ibridazione competitiva, i due canali nel nostro microarray Agilent hanno dimostrato di essere relativamente indipendente come saturazione in uno dei canali è molto raro. Di tutti i profili INPUT-chip che abbiamo estratto, presentiamo solo l’analisi di otto profili rappresentativi qui (due da ciascuno dei quattro punti di tempo di sviluppo) in quanto la maggior parte dei profili INPUT-chip sono molto simili (Additional file 2: Figura S1). Gli otto profili INPUT-chip sono stati poi confrontati con i nove profili INPUT-seq raccolti in questo studio (Additional file 1: Tabella S3).

Una delle osservazioni più sorprendenti è che i profili INPUT-chip e INPUT-seq sembrano essere sostanzialmente diversi, anche se lo stesso materiale DNA di input è stato utilizzato per l’ibridazione microarray e sequenziamento (Figura 1). La grandezza relativa e la posizione dei picchi sembrano essere coerenti attraverso i profili INPUT-chip da più esperimenti. Tuttavia, i modelli nei nove profili INPUT-seq sembrano essere più variabili. Possiamo identificare visivamente molte regioni che hanno incoerente arricchimento del segnale attraverso più profili INPUT-seq (evidenziato nella Figura 1a). Un’analisi di clustering è stata eseguita per quantificare questa osservazione. Abbiamo trovato che tutti gli otto profili INPUT-chip raggruppati strettamente a vicenda (Figura 1b). Questo risultato dimostra che la distribuzione del DNA di fondo misurata da microarray e sequenziamento ad alta velocità è diverso. Tutti i profili INPUT-chip e sette dei nove profili INPUT-seq correlati positivamente con il contenuto genomico GC a livello genomico (Figura 1b), così come intorno ai siti di inizio trascrizione (TSS) e siti di fine trascrizione (TES) (Figura 1c). La forza della correlazione con GC è altamente coerente tra i profili INPUT-chip, ma altamente variabile tra i profili INPUT-seq (Figura 1b-c e file aggiuntivo 2: Figura S2). In particolare, i profili INPUT-seq ottenuti a E-16-20 h (E16) e E-20-24 h (E20) non correlano con il contenuto di GC.

Si nota anche che INPUT-seq con maggiore profondità di sequenziamento (>4 milioni di letture mappate) tendono a raggrupparsi più strettamente di quelli con profondità di sequenziamento inferiore, suggerendo che ci può essere una relazione tra profondità di sequenziamento e variabilità del DNA di ingresso. Per testare questa ipotesi, abbiamo generato 11 ulteriori profili INPUT-seq sottocampionando le letture di sequenziamento dal campione di DNA di ingresso più profondamente sequenziato (AdultMale; AM) a diverse proporzioni di campionamento (Figura 1d e file aggiuntivo 2: Figura S3). Come previsto, i profili con maggiore profondità di sequenziamento tendono a raggrupparsi più fortemente insieme, e la loro correlazione con la variazione del contenuto di GC è più coerente. Tuttavia, la correlazione del contenuto di GC diventa molto più debole solo ad una profondità di sequenziamento molto bassa (<2 milioni di letture; Figura 1d). Ciò indica che la bassa profondità di sequenziamento non è l’unico fattore che influenza la qualità INPUT-seq. Inoltre, alcuni INPUT-seq con profondità di sequenziamento relativamente bassa (E0 e AF, <4 milioni di letture) può dare profili di DNA di input coerenti. Questo implica che INPUT-seq variabilità può anche essere attribuito ad altri fattori sperimentali. Anche se ulteriori studi sono necessari per sezionare l’intera gamma di fattori sperimentali che influenzano la variabilità delle librerie di DNA di ingresso, potrebbe essere influenzato da variazioni nella preparazione del campione (ad esempio, la preparazione cromatina diversa e sonicazione), run-to-run variazione del sequenziatore, sequencer-to-sequencer variazione per lo stesso modello, e una serie di altre variabili negli esperimenti. L’alta variabilità tra i profili INPUT-seq è davvero un problema critico, poiché la grande variabilità contribuisce all’instabilità della stima della densità in un profilo ChIP-seq, che influenzerà l’analisi dei dati a valle. Come verrà mostrato nelle sezioni successive di questo documento, un INPUT-seq con correlazione insolitamente debole con il contenuto di GC può avere un impatto sulla costruzione di profili medi in posizioni genomiche importanti. E ‘quindi indispensabile per sequenziare il DNA di ingresso a profondità sufficiente e per accertare che il profilo ottenuto è coerente con quelli da esperimenti simili.

Copertura genomica è un’altra considerazione chiave quando si sceglie tra ChIP-chip e ChIP-seq. La copertura genomica di ChIP-chip è limitata dal disegno della sonda microarray, e la copertura di ChIP-seq dipende dalla profondità di sequenziamento. La copertura genomica raggiunta dal nostro microarray Agilent è di circa il 70%. Utilizzando i dati INPUT-seq sotto-campionati, dimostriamo che INPUT-seq fornisce generalmente una copertura genomica più alta a profondità di sequenziamento bassa come un milione di letture. Questa tendenza costruita dai dati sotto-campionati in modo casuale conferma la copertura genomica osservata degli altri otto set di dati INPUT-seq reali (Figura 1e).

Confronto delle caratteristiche del profilo

Abbiamo poi confrontato le caratteristiche dei profili ChIP-chip e ChIP-seq. Per confrontare i profili generati dalle due tecnologie, abbiamo diviso il genoma in 1 kb non sovrapposti bin e definito il livello di arricchimento in ogni bin come la media del rapporto log del canale IP sul canale di ingresso (vedi la sezione Metodi per i dettagli). Ci riferiamo alla distribuzione del segnale di un profilo ChIP come la sua distribuzione dei valori di arricchimento di tutti i bin. In primo luogo, abbiamo mirato a caratterizzare il rapporto medio segnale-rumore per i profili generati dalle due tecnologie. Abbiamo usato l’asimmetria (troncata) del profilo di densità del segnale dopo aver rimosso i segnali dal 5% più alto e più basso della distribuzione come misura del rapporto segnale-rumore di un profilo. L’asimmetria è una misura dell’asimmetria di una distribuzione e un’asimmetria positiva indica che la coda sul lato destro è più lunga, il che implica un buon rapporto segnale-rumore. In quasi tutti i casi, un profilo ChIP-seq ha una skewness più alta rispetto al suo corrispondente profilo ChIP-chip della stessa condizione biologica (Figura 2 e Additional file 1: Tabella S4). Notiamo che la differenza di skewness dipende dal fattore IP che potrebbe essere dovuto alla diversa qualità degli anticorpi e alla prevalenza di eventi di modifica o di legame degli istoni. La stessa conclusione può essere tratta anche se è stata utilizzata una dimensione diversa bin (file aggiuntivo 2: Figura S4). I nostri risultati hanno confermato l’osservazione generale che ChIP-seq di solito produce un profilo di segnale più distintivo di ChIP-chip.

In seguito, abbiamo caratterizzato le regioni di arricchimento all’interno di ogni profilo ChIP. Per eseguire un confronto equo, vorremmo utilizzare un algoritmo che esegue il peak calling sui dati ChIP-seq e ChIP-chip utilizzando gli stessi criteri. Attualmente, molti algoritmi di chiamata di picco comunemente usati sono specificamente progettati per analizzare i dati ChIP-chip o ChIP-seq, ma non entrambi. Per superare questa limitazione, abbiamo identificato i picchi da entrambi i profili ChIP-chip e ChIP-seq utilizzando lo stesso genoma-scansione euristica (vedi la sezione Metodi). I nostri risultati indicano che possiamo quasi sempre scoprire un maggior numero di picchi e picchi più stretti utilizzando i dati generati da ChIP-seq rispetto a ChIP-chip quando si analizza lo stesso campione biologico, e questa conclusione è coerente indipendentemente dalla rigorosità dei criteri di identificazione utilizzati (Figura 2 e Additional file 2: Figura S5). In pratica, possiamo probabilmente identificare un numero ancora maggiore di picchi stretti nei dati ChIP-seq se facciamo esplicitamente uso di informazioni specifiche del filamento all’interno della procedura di chiamata di picco (oltre a spostare solo ogni lettura verso la sua estremità 5 ‘di un numero costante di coppia di basi), quindi l’analisi attuale fornisce un limite inferiore sull’efficacia di ChIP-seq rispetto a ChIP-chip. Nel complesso, i nostri risultati dimostrano che ChIP-seq fornisce una maggiore risoluzione spaziale e il rapporto segnale-rumore.

Riproducibilità del segnale a livello genomico all’interno e tra le tecnologie

Inoltre, abbiamo stimato la riproducibilità tra profili ChIP-chip e/o ChIP-seq a livello genomico (1 kb bins). Per evitare distorsioni dovute a differenze nella copertura genomica e nella mappatura della sequenza (Figura 1e), abbiamo escluso le regioni genomiche che non contengono sonde microarray e le regioni con una variabilità insolitamente alta attraverso più profili INPUT-seq. Il coefficiente di correlazione di Pearson, r, è stato utilizzato come misura di correlazione, poiché è più sensibile del coefficiente di correlazione di Spearman per confrontare la coda di due distribuzioni di segnale, che è particolarmente importante nell’analisi dei profili di segnale di arricchimento ChIP. La correlazione tra coppie di repliche ChIP-chip e tra coppie di repliche ChIP-seq è generalmente alta (mediana r = 0,85 e 0,82, rispettivamente), indicando che entrambe le tecnologie possono produrre risultati riproducibili. Come previsto, la correlazione cross-platform tra coppie replicate di profili ChIP-chip e ChIP-seq sono più modeste (mediana r = 0,41; Additional file 1: Tabella S5). Conclusioni simili possono essere raggiunte anche se usiamo diverse dimensioni bin per calcolare la correlazione inter-profilo (Additional file 2: Figura S6). Un grafico di dispersione rappresentativo che confronta ogni coppia di tecnologie è mostrato nella Figura 3b-d. Osserviamo anche una correlazione positiva tra la skewness e la riproducibilità inter-profilo (Additional file 2: Figura S7), suggerendo anticorpi più sensibili possono produrre profili più coerenti tra le due tecnologie.

Costruzione del profilo del segnale medio a TSS e TES

Costruzione di profili di segnale ChIP media intorno a importanti caratteristiche genomiche come TSS e TES è un modo comune per visualizzare l’arricchimento del segnale intorno a queste caratteristiche. Pertanto, abbiamo studiato la riproducibilità della media TSS e profili TES (2 kb a monte e 2 kb a valle) per ogni coppia di profili ChIP replicati (Additional file 2: Figura S8). I profili medi della maggior parte delle coppie replicate sono altamente coerenti. Tuttavia, ci sono alcune coppie che sono significativamente diversi, in particolare i profili di H3K27Me3 e H3K9Me3 sia in fase E-16-20 h e E-20-24 h (Additional file 2: Figure S8c e S8g). Senza convalida esterna, è impossibile determinare se i profili di segnale medio generato da ChIP-chip o ChIP-seq sono più accurate. Tuttavia, due linee di prova ci hanno portato a credere che i profili di segnale medio da ChIP-chip erano più probabilità di essere accurate. In primo luogo, tutte e tre le repliche ChIP-chip in questi punti temporali avevano profili medi molto coerenti. In secondo luogo, i profili di segnale medio ChIP-seq in queste condizioni biologiche assomigliavano alla tendenza della variazione del contenuto di GC a TSS e TES (Figura 1c). Le correlazioni insolitamente basse tra i contenuti GC e i profili INPUT-seq di E-16-20 h e E-20-24 h (Figura 1b e file aggiuntivo 2: Figura S2b) ci ha spinto a ipotizzare che la discrepanza osservata era dovuta a una rappresentazione errata della variazione del contenuto GC dai rispettivi profili INPUT-seq. Sia H3K27Me3 e H3K9Me3 sono marchi repressivi che di solito sono impoveriti a TSSs e TESs così qualsiasi variazione nella sottrazione di fondo è probabilmente molto più pronunciato rispetto ad altri marchi istone che hanno forte arricchimento del segnale a queste caratteristiche genomiche. Per testare la nostra ipotesi, abbiamo sostituito lo sfondo INPUT-seq corrispondente con l’INPUT-seq dal campione AdultFemale, poiché ha la più alta correlazione con la variazione del contenuto di GC. Dopo la sostituzione, i profili medi dei segnali generati da ChIP-seq e ChIP-chip in questi due stadi di sviluppo sono d’accordo (Figura 4 e file aggiuntivo 2: Figura S9). Questo risultato è sorprendente in quanto dimostra che utilizzando diversi INPUT-seq come controllo negativo dello stesso profilo ChIP-seq può portare a un’interpretazione sostanzialmente diversa dei dati.

Effetto dell’uso di diversi profili di input nella normalizzazione dei dati ChIP-seq

Avendo osservato l’impatto di INPUT-seq nella costruzione di profili TSS e TES medi, ci siamo chiesti se l’uso di diversi profili INPUT-seq per la normalizzazione dello sfondo influisce significativamente sui risultati di ChIP-seq. Abbiamo usato SPP per chiamare i picchi per 10 dei nostri campioni ChIP-seq (CBP, H3K9Ac, H3K9Me3, H3K27Ac, H3K27Me3 a E16-20 h e E20-24 h) dove ogni profilo ChIP è stato normalizzato contro quattro diversi INPUT-seq come sfondo (l’ingresso dal punto di tempo corrispondente, AdultFemale, AdultMale, e E-4-8 h). Questi profili INPUT-seq sono stati scelti perché hanno diversa profondità di sequenziamento e la correlazione con il contenuto di GC (Figura 1b). Un confronto tra il numero di picchi e la larghezza mediana del picco è mostrato nella Figura 5. Abbiamo osservato una grande differenza nel numero di picchi di essere chiamato per qualsiasi campione ChIP-seq quando diversi INPUT-seq è stato utilizzato come sfondo. Nel caso estremo (E-16-24 h, H3K9Me3 ChIP), il numero di picchi può cambiare da zero a quasi 40.000 ad un FDR del 5% (Figura 5a). In generale, più picchi statisticamente significativi (FDR < 0,05) sono stati rilevati quando si normalizza contro un campione di DNA di input profondamente sequenziato (AdultMale e E-4-8 h in questo esperimento), anche se la grandezza assoluta della differenza varia tra i set di dati ChIP. La differenza nel numero di picchi indica probabilmente una differenza nel potere di rilevamento. Per ogni campione ChIP, abbiamo calcolato la proporzione di sovrapposizione tra ogni coppia di set di picchi generati da quattro diversi background DNA di ingresso (cioè, sei confronti per campione ChIP). Abbiamo trovato che la proporzione media di sovrapposizione rispetto al set di picchi più piccolo è di circa il 95%, indicando che le differenze nel numero di picco rilevato è probabilmente dovuto al diverso potere di chiamare i picchi più deboli. Abbiamo osservato che i picchi forti (cioè, quelli con basso rilevamento FDR) sono stati più probabilmente rilevato in diversi set di picchi (vedi file aggiuntivo 2: Figura S10 per un esempio). La larghezza mediana dei picchi rilevati è anche influenzato utilizzando diversi INPUT-seq come sfondo (Figura 5b). Questa analisi ha dimostrato che la normalizzazione utilizzando diversi INPUT-seq può avere un impatto significativo, e sottovalutato, sulla chiamata di picco.

Valutazione della variazione dovuta all’uso di diversi peak callers

Un’altra importante fonte di variazione nell’analisi dei profili ChIP-chip e ChIP-seq deriva dall’uso di diversi algoritmi di analisi. Un gran numero di strumenti di analisi ChIP-chip e ChIP-seq disponibili al pubblico sono stati sviluppati fino ad oggi, e tutti utilizzano diversi metodi per il tag shifting, la normalizzazione del profilo, lo smoothing, l’identificazione dei picchi e il calcolo del false discovery rate. Non è quindi troppo sorprendente trovare che i diversi chiamanti di picco possono generare risultati molto diversi in termini di identificazione del sito di legame, in particolare quando si tratta di picchi con segnali deboli. Utilizzando il nostro compendio di dataset ChIP-chip e ChIP-seq, abbiamo potuto valutare quanta variazione nell’identificazione dei picchi può essere attribuita all’uso di diverse tecnologie di profilazione e all’uso di diversi peak caller. In questo studio, abbiamo analizzato i nostri profili ChIP-chip utilizzando due peak caller: MA2C e Splitter e abbiamo analizzato i nostri profili ChIP-seq utilizzando altri due peak callers: MACS e SPP (vedi file aggiuntivo 1: Tabella S8). Questi chiamatori di picco sono stati scelti perché sono ampiamente utilizzati, pubblicamente disponibili, e generalmente mostrano buone prestazioni in precedenti studi comparativi. Abbiamo calcolato la sovrapposizione dei primi 1.000 picchi di quattro dei fattori (CBP, H3K4Me1, H3K4Me3, e H3K27Me3) attraverso più fasi di sviluppo. I quattro fattori IP sono stati scelti in quanto erano profili rappresentativi contenenti picchi ampi (CBP e H3K27Me3) e picchi stretti (H3K4Me1 e H3K4Me3). Qui, presentiamo solo i risultati del confronto dei primi 1.000 picchi, poiché questo è un numero biologicamente ragionevole di siti di arricchimento ad alta confidenza in questi profili. La conclusione generale di questa analisi è robusto contro una varietà di picco chiamando soglie (Additional file 2: Figura S11). La concordanza tra due set di picchi è stata misurata dalla proporzione media di picchi sovrapposti. Come mostrato nella Figura 6, i confronti basati su profili di H3K4Me1 e H3K4Me3 prodotto risultati attesi, in cui la concordanza intra-piattaforma è superiore a cross-piattaforma concordanza (cioè, set di picco generato da due chiamanti di picco sul profilo stesso sono più concordanti di set di picco generato da due chiamanti di picco su due profili). Tuttavia, la concordanza intra-piattaforma può essere bassa come la concordanza inter-piattaforma quando si analizzano i profili di H3K27Me3 e CBP, il che implica che la variazione negli algoritmi di chiamata di picco può essere grande come l’uso di diverse tecnologie di profiling per alcuni fattori IP. L’osservazione che gli attuali algoritmi di chiamata di picco producono risultati meno concordanti per i profili ChIP con ampi domini (CBP e H3K27Me3) rispetto a quelli con picchi netti (H3K4Me1 e H3K4Me3) può suggerire che sono meno coerenti nell’identificare ampie regioni di arricchimento, che può essere un argomento interessante per ulteriori indagini.