ChIP-chip versus ChIP-seq: Tanulságok a kísérleti tervezéshez és az adatelemzéshez

A bemeneti DNS-profilok elemzése

A szekvenáláson alapuló és a microarray-alapú ChIP-adatok közötti technológiai különbségek megértéséhez először a microarray (INPUT-chip) és a nagy áteresztőképességű szekvenálás (INPUT-seq) által létrehozott, keresztkötésű és szonikált DNS-fragmentumok (bemeneti DNS) profiljait elemeztük. Mivel a bemeneti DNS-profilnak függetlennek kell lennie a ChIP-hez használt antitesttől, ez az összehasonlítás betekintést nyújthat a két profilkészítési technológia közötti specifikus különbségekbe. Az INPUT-chip adatokat a kétcsatornás microarray-adataink háttércsatornájából nyertük. Míg ez a microarray platform kompetitív hibridizációt használ, a mi Agilent microarray-nk két csatornája viszonylag függetlennek bizonyult, mivel a telítődés bármelyik csatornában nagyon ritka . Az összes kivont INPUT-chip profil közül itt csak nyolc reprezentatív profil elemzését mutatjuk be (kettőt mind a négy fejlődési időpontból), mivel az INPUT-chip profilok többsége nagyon hasonló (Additional file 2: S1 ábra). A nyolc INPUT-chip-profilt ezután összehasonlítottuk az ebben a vizsgálatban gyűjtött kilenc INPUT-seq-profillal (Additional file 1: S3 táblázat).

Az egyik legszembetűnőbb megfigyelés, hogy az INPUT-chip és az INPUT-seq profilok lényegesen különbözőnek tűnnek, annak ellenére, hogy a microarray-hibridizációhoz és a szekvenáláshoz ugyanazt az input DNS-anyagot használtuk (1. ábra). A csúcsok relatív nagysága és elhelyezkedése konzisztensnek tűnik a több kísérletből származó INPUT-chip profilokban. A kilenc INPUT-seq-profilban azonban a minták sokkal változatosabbnak tűnnek. Vizuálisan számos olyan régiót azonosíthatunk, amelyek több INPUT-seq-profilban is következetlen jelgazdagodást mutatnak (az 1a. ábrán kiemelve). Ennek a megfigyelésnek a számszerűsítésére klaszterezési elemzést végeztünk. Azt találtuk, hogy mind a nyolc INPUT-chip profil szorosan egymáshoz klasztereződött (1b. ábra). Ez az eredmény azt mutatja, hogy a microarray és a nagy áteresztőképességű szekvenálás alapján mért háttér-DNS-eloszlás eltérő. Az összes INPUT-chip és a kilenc INPUT-seq profilból hét pozitívan korrelált a genom GC-tartalmával genomszintű szinten (1b. ábra), valamint a transzkripciós starthelyek (TSS) és a transzkripciós véghelyek (TES) körül (1c. ábra). A GC-vel való korreláció erőssége nagyon konzisztens az INPUT-chip profilok között, de nagyon változó az INPUT-seq profilok között (1b-c. ábra és Additional file 2: S2. ábra). Figyelemre méltó, hogy az E-16-20 h-ban (E16) és az E-20-24 h-ban (E20) kapott INPUT-seq profilok nem korrelálnak a GC-tartalommal.

Megjegyezzük azt is, hogy a nagyobb szekvenálási mélységű (>4 millió leképezett olvasat) INPUT-seq-ek általában szorosabban klasztereződnek, mint az alacsonyabb szekvenálási mélységűek, ami arra utal, hogy kapcsolat lehet a szekvenálási mélység és az input DNS variabilitása között. E hipotézis tesztelésére 11 további INPUT-seq-profilt hoztunk létre a legmélyebben szekvenált bemeneti DNS-mintából (AdultMale; AM) származó szekvenálási olvasatok különböző mintavételi arányú almintavételezésével (1d. ábra és Additional file 2: S3 ábra). A várakozásoknak megfelelően a nagyobb szekvenálási mélységű profilok általában erősebben klasztereződnek, és korrelációjuk a GC-tartalom változásával konzisztensebb. A GC-tartalom korreláció azonban csak nagyon alacsony szekvenálási mélységnél (<2 millió leolvasás; 1d. ábra) válik sokkal gyengébbé. Ez azt jelzi, hogy az alacsony szekvenálási mélység nem az egyetlen tényező, amely befolyásolja az INPUT-seq minőségét. Sőt, néhány viszonylag alacsony szekvenálási mélységű INPUT-seq (E0 és AF, <4 millió olvasat) is adhat konzisztens bemeneti DNS-profilokat. Ez azt jelenti, hogy az INPUT-seq variabilitása más kísérleti tényezőknek is tulajdonítható. Bár további vizsgálatokra van szükség a bemeneti DNS-könyvtárak variabilitását befolyásoló kísérleti tényezők teljes körének feltárásához, azt befolyásolhatják a minta előkészítésének eltérései (pl. különböző kromatin-előkészítés és szonikálás), a szekvenáló futásonkénti eltérései, ugyanazon modell esetében a szekvenáló és szekvenáló közötti eltérések, valamint a kísérletek számos más változója. Az INPUT-seq-profilok közötti nagyfokú variabilitás valóban kritikus kérdés, mivel a nagyfokú variabilitás hozzájárul a ChIP-seq-profil sűrűségbecslésének instabilitásához, ami befolyásolja a későbbi adatelemzést. Amint azt e cikk későbbi szakaszaiban bemutatjuk, a GC-tartalommal szokatlanul gyenge korrelációt mutató INPUT-seq befolyásolhatja az átlagos profilok felépítését fontos genomi helyeken. Ezért elengedhetetlen a bemeneti DNS megfelelő mélységű szekvenálása, és annak biztosítása, hogy a kapott profil összhangban legyen a hasonló kísérletekből származó profilokkal.

A genomiális lefedettség egy másik kulcsfontosságú szempont a ChIP-chip és a ChIP-seq közötti választás során. A ChIP-chip genomi lefedettségét a microarray szondák kialakítása korlátozza, a ChIP-seq lefedettsége pedig a szekvenálási mélységtől függ. Az Agilent microarray-vel elért genomiális lefedettség körülbelül 70%. Az almintavételezett INPUT-seq adatok felhasználásával megmutatjuk, hogy az INPUT-seq általában nagyobb genomi lefedettséget biztosít már egymillió leolvasás mélységű szekvenálási mélységnél is. Ez a véletlenszerűen almintavételezett adatokból konstruált tendencia megerősíti a másik nyolc valódi INPUT-seq-adatkészlet megfigyelt genomi lefedettségét (1e. ábra).

Profiljellemzők összehasonlítása

Ezután összehasonlítottuk a ChIP-chip és a ChIP-seq profilok jellemzőit. A két technológia által generált profilok összehasonlításához a genomot 1 kb-os, nem átfedő binsekre osztottuk, és az egyes binekben a dúsulási szintet az IP-csatorna és a bemeneti csatorna logaritmikus arányának átlagaként határoztuk meg (a részleteket lásd a Módszerek részben). A ChIP-profil jeleloszlására úgy hivatkozunk, mint az összes bins dúsulási értékének eloszlására. Először a két technológia által generált profilok átlagos jel-zaj arányának jellemzésére törekedtünk. Egy profil jel-zaj arányának mérésére a jelsűrűségprofil (csonka) ferdeségét használtuk, miután eltávolítottuk az eloszlás legmagasabb és legalacsonyabb 5%-ából származó jeleket. A ferdeség az eloszlás aszimmetriájának mérőszáma, és a pozitív ferdeség azt jelzi, hogy a jobb oldali farok hosszabb, ami jó jel-zaj arányra utal. Szinte minden esetben a ChIP-seq profil nagyobb ferdeséggel rendelkezik, mint a megfelelő ChIP-chip profil ugyanabban a biológiai állapotban (2. ábra és Additional file 1: Table S4). Megjegyezzük, hogy a ferdeség különbsége az IP-tényezőtől függ, ami az eltérő antitestminőségnek és a hisztonmódosulások vagy kötődési események előfordulási gyakoriságának tudható be. Ugyanez a következtetés vonható le akkor is, ha eltérő bin méreteket használtunk (Additional file 2: Figure S4). Eredményeink megerősítették azt az általános megfigyelést, hogy a ChIP-seq általában markánsabb jelprofilt eredményez, mint a ChIP-chip.

A következőkben jellemeztük a gazdagodási régiókat az egyes ChIP-profilokon belül. A tisztességes összehasonlítás érdekében olyan algoritmust szeretnénk használni, amely a ChIP-seq és a ChIP-chip adatokon azonos kritériumok alapján végzi a csúcsmeghívást. Jelenleg sok általánosan használt csúcshívó algoritmus kifejezetten ChIP-chip vagy ChIP-seq adatok elemzésére készült, de nem mindkettőre. Ennek a korlátozásnak a kiküszöbölése érdekében a ChIP-chip és a ChIP-seq profilokból azonos genompásztázó heurisztika segítségével azonosítottuk a csúcsokat (lásd a Módszerek fejezetet). Eredményeink azt mutatják, hogy ugyanazon biológiai minta elemzésekor szinte mindig nagyobb számú csúcsot és keskenyebb csúcsokat fedezhetünk fel a ChIP-seq-ről generált adatok felhasználásával, mint a ChIP-chip-ről, és ez a következtetés az alkalmazott azonosítási kritériumok szigorától függetlenül konzisztens (2. ábra és Additional file 2: Figure S5). A gyakorlatban valószínűleg még nagyobb számú keskeny csúcsot tudunk azonosítani a ChIP-seq-adatokban, ha a csúcsmeghatározási eljáráson belül kifejezetten felhasználjuk a szálspecifikus információkat (amellett, hogy minden olvasatot csak az 5′ vége felé tolunk el egy állandó bázispárszámmal), így a jelenlegi elemzés egy alsó korlátot ad a ChIP-seq hatékonyságára a ChIP-chiphez képest. Összességében eredményeink azt mutatják, hogy a ChIP-seq nagyobb térbeli felbontást és jel-zaj arányt biztosít.

Genomszintű jel reprodukálhatóság a technológiákon belül és a technológiák között

A továbbiakban a ChIP-chip és/vagy ChIP-seq profilok közötti reprodukálhatóságot genomszintű szinten (1 kb bins) becsültük meg. A genomi lefedettség és a szekvenciatérképezés különbségeiből eredő torzítások elkerülése érdekében (1e. ábra) kizártuk azokat a genomi régiókat, amelyek nem tartalmaznak mikroarray-szondákat, valamint azokat a régiókat, amelyek szokatlanul nagy variabilitást mutatnak több INPUT-seq-profilon keresztül. A Pearson-féle korrelációs együtthatót (r) használtuk a korreláció mérésére, mivel ez érzékenyebb, mint a Spearman-féle korrelációs együttható két jeleloszlás farokrészének összehasonlítására, ami különösen fontos a ChIP-dúsítási jelprofilok elemzésénél. A ChIP-chip-replikátumpárok és a ChIP-seq-replikátumpárok közötti korreláció általában magas (medián r = 0,85, illetve 0,82), ami azt jelzi, hogy mindkét technológia reprodukálható eredményeket produkál. A várakozásoknak megfelelően a ChIP-chip és a ChIP-seq profilok replikátumpárjai közötti keresztplatformos korreláció szerényebb (medián r = 0,41; Additional file 1: Table S5). Hasonló következtetésekre juthatunk akkor is, ha a profilok közötti korreláció kiszámításához különböző bin méreteket használunk (Additional file 2: Figure S6). Az egyes technológiai párokat összehasonlító reprezentatív szórásdiagramot a 3b-d. ábra mutatja. Pozitív korrelációt figyelhetünk meg a ferdeség és a profilok közötti reprodukálhatóság között is (Additional file 2: Figure S7), ami arra utal, hogy az érzékenyebb antitestek konzisztensebb profilokat eredményezhetnek a két technológia között.

A TSS és TES körüli átlagos jelprofil konstruálása

A fontos genomi jellemzők, például a TSS és a TES körüli átlagos ChIP-jelprofilok konstruálása gyakori módja az e jellemzők körüli jelgazdagodás megjelenítésének. Ezért megvizsgáltuk az átlagos TSS és TES profilok (2 kb felfelé és 2 kb lefelé) reprodukálhatóságát minden pár replikált ChIP-profil esetében (Additional file 2: S8 ábra). A legtöbb replikátumpár átlagos profiljai nagymértékben konzisztensek. Van azonban néhány pár, amely jelentősen eltér, különösen a H3K27Me3 és a H3K9Me3 profiljai az E-16-20 h és az E-20-24 h szakaszban (Additional file 2: S8c és S8g ábra). Külső validálás nélkül nem lehet eldönteni, hogy a ChIP-chip vagy a ChIP-seq által generált átlagos jelprofilok pontosabbak-e. Mindazonáltal két bizonyíték arra engedett következtetni, hogy a ChIP-chipből származó átlagos jelprofilok nagyobb valószínűséggel pontosabbak. Először is, mindhárom ChIP-chip-replikátum ezekben az időpontokban nagyon konzisztens átlagos profilokat mutatott. Másodszor, a ChIP-seq átlagos jelprofiljai ezeknél a biológiai körülményeknél hasonlítottak a GC-tartalom változásának trendjéhez a TSS és a TES esetében (1c. ábra). A GC-tartalom és az E-16-20 h és E-20-24 h INPUT-seq-profilok közötti szokatlanul alacsony korrelációk (1b. ábra és Additional file 2: S2b ábra) arra késztettek minket, hogy feltételezzük, hogy a megfigyelt eltérés a GC-tartalom változásának az adott INPUT-seq-profilok általi téves ábrázolásából ered. Mind a H3K27Me3, mind a H3K9Me3 represszív jelölések, amelyek általában a TSS-eknél és a TES-eknél kimerülnek, így a háttérkivonás során tapasztalt bármilyen eltérés valószínűleg sokkal hangsúlyosabb, mint más olyan hiszton-jelölések esetében, amelyek erős jelgazdagodással rendelkeznek ezeknél a genomiális jellemzőknél. Hipotézisünk teszteléséhez a megfelelő INPUT-seq hátteret az AdultFemale mintából származó INPUT-seq-re cseréltük, mivel ez rendelkezik a legmagasabb korrelációval a GC-tartalom változásával. A csere után a ChIP-seq és a ChIP-chip által generált átlagos jelprofilok ebben a két fejlődési szakaszban megegyeznek (4. ábra és Additional file 2: S9 ábra). Ez az eredmény azért figyelemre méltó, mert azt mutatja, hogy ugyanazon ChIP-seq-profil különböző INPUT-seq negatív kontrolljaként való használata lényegesen eltérő adatértelmezéshez vezethet.

A különböző bemeneti profilok használatának hatása a ChIP-seq-adatok normalizálásában

Miután megfigyeltük az INPUT-seq hatását az átlagos TSS- és TES-profilok felépítésében, megkérdeztük, hogy a különböző INPUT-seq-profilok használata a háttér normalizálásához jelentősen befolyásolja-e a ChIP-seq-csúcshívási eredményeket. SPP-t használtunk 10 ChIP-seq-mintánk (CBP, H3K9Ac, H3K9Me3, H3K27Ac, H3K27Me3 az E16-20 h és E20-24 h időpontban) csúcsok hívására, ahol minden egyes ChIP-profilt háttérként négy különböző INPUT-seq-hez képest normalizáltunk (a megfelelő időpontból származó input, AdultFemale, AdultMale és E-4-8 h). Ezeket az INPUT-seq profilokat azért választottuk, mert különböző szekvenálási mélységgel és GC-tartalommal korrelálnak (1b. ábra). A csúcsok számának és a csúcsok medián szélességének összehasonlítása az 5. ábrán látható. Nagy különbséget figyeltünk meg a hívott csúcsok számában bármely ChIP-seq minta esetében, amikor különböző INPUT-seq-et használtunk háttérként. A szélsőséges esetben (E-16-24 h, H3K9Me3 ChIP) a csúcsok száma 5%-os FDR mellett nullától közel 40 000-ig változhat (5a. ábra). Általában több statisztikailag szignifikáns csúcsot (FDR < 0,05) észleltünk, amikor mélyen szekvenált bemeneti DNS-mintával (AdultMale és E-4-8 h ebben a kísérletben) szemben normalizáltunk, bár a különbség abszolút nagysága ChIP-adatkészletek között változik. A csúcsok számában mutatkozó különbség valószínűleg a detektálási teljesítmény különbségét jelzi. Minden egyes ChIP-mintára kiszámítottuk az átfedés arányát a négy különböző bemeneti DNS-háttérrel létrehozott csúcskészletek minden egyes párja között (azaz ChIP-mintánként hat összehasonlítás). Azt találtuk, hogy az átfedés átlagos aránya a kisebb csúcskészlet tekintetében körülbelül 95%, ami azt jelzi, hogy a detektált csúcsok számában mutatkozó különbségek valószínűleg a gyengébb csúcsok hívásának eltérő teljesítményéből adódnak. Megfigyeltük, hogy az erős csúcsokat (azaz az alacsony detektálási FDR-rel rendelkező csúcsokat) nagyobb valószínűséggel detektálták a különböző csúcskészletekben (lásd a 2. kiegészítő fájl: S10. ábra példáját). A detektált csúcsok medián szélességét az is befolyásolja, ha különböző INPUT-seq-et használunk háttérként (5b. ábra). Ez az elemzés azt mutatta, hogy a különböző INPUT-seq-et használó normalizálásnak jelentős és alulértékelt hatása lehet a csúcsmeghatározásra.

A különböző csúcshívók használatából eredő eltérések értékelése

A ChIP-chip és ChIP-seq profilok elemzésének másik fontos eltérési forrása a különböző elemzési algoritmusok használatából ered. A mai napig számos nyilvánosan elérhető ChIP-chip és ChIP-seq elemző eszközt fejlesztettek ki, és ezek mindegyike különböző módszereket alkalmaz a címkék eltolására, a profil normalizálására, a simításra, a csúcsok azonosítására és a hamis felfedezési arány kiszámítására. Ezért nem túl meglepő, hogy a különböző csúcshívók meglehetősen eltérő eredményeket hozhatnak a kötőhelyek azonosítása tekintetében, különösen a gyenge jelekkel rendelkező csúcsok esetében. A ChIP-chip és ChIP-seq adatkészletekből összeállított összeállításunk segítségével felmérhettük, hogy a csúcsok azonosításában mutatkozó eltérések mekkora része vezethető vissza a különböző profilalkotási technológiák és a különböző csúcshívók használatára. Ebben a tanulmányban a ChIP-chip-profiljainkat két peak caller segítségével elemeztük: MA2C és Splitter, a ChIP-seq-profiljainkat pedig két másik csúcshívó segítségével elemeztük: MACS és SPP (lásd 1. kiegészítő fájl: S8. táblázat). Ezeket a csúcshívókat azért választottuk, mert széles körben használtak, nyilvánosan elérhetőek, és általában jó teljesítményt mutattak korábbi összehasonlító vizsgálatokban . Négy faktor (CBP, H3K4Me1, H3K4Me3, H3K4Me3 és H3K27Me3) felső 1000 csúcsának átfedését számoltuk ki több fejlődési szakaszban. A négy IP-faktort azért választottuk, mert reprezentatív profilok voltak, amelyek széles csúcsokat (CBP és H3K27Me3) és keskeny csúcsokat (H3K4Me1 és H3K4Me3) tartalmaztak. Itt csak a legfelső 1000 csúcs összehasonlításának eredményeit mutatjuk be, mivel ez biológiailag ésszerű számú, magas konfidenciájú dúsulási helyek száma ezekben a profilokban. Ennek az elemzésnek az általános következtetése robusztus a különböző csúcshívási küszöbértékekkel szemben (Additional file 2: S11 ábra). A két csúcskészlet közötti összhangot az átfedő csúcsok átlagos arányával mértük. Amint az a 6. ábrán látható, a H3K4Me1 és H3K4Me3 profilokon alapuló összehasonlítások a várt eredményeket hozták, amelyek szerint a platformon belüli konkordancia magasabb, mint a platformok közötti konkordancia (azaz a két csúcshívó által ugyanazon profilon létrehozott csúcskészletek jobban egyeznek, mint a két csúcshívó által két profilon létrehozott csúcskészletek). A H3K27Me3 és a CBP profiljainak elemzésekor azonban a platformon belüli konkordancia olyan alacsony lehet, mint a platformok közötti konkordancia, ami arra utal, hogy a csúcshívó algoritmusok eltérései ugyanolyan nagyok lehetnek, mint a különböző profilalkotási technológiák használata egyes IP-faktorok esetében. Az a megfigyelés, hogy a jelenlegi csúcshívó algoritmusok kevésbé egybehangzó eredményeket produkálnak a széles doméneket (CBP és H3K27Me3) tartalmazó ChIP-profilok esetében, mint az éles csúcsokat (H3K4Me1 és H3K4Me3) tartalmazó profilok esetében, arra utalhat, hogy kevésbé következetesek a széles dúsulási régiók azonosításában, ami további vizsgálatok érdekes témája lehet.