Experimental procedure
Patiens xenografting and cell culture
A humán T-ALL sejtvonalak közé tartozik a CUTLL1 (ajándék Adolfo Ferrando, Columbia Egyetem) és a JURKAT (American Type Culture Collection (ATCC), Manassas, VA, #CCL-119) . A sejteket RPMI1640 médiumban tenyésztettük, l-glutaminnal és 25 mM HEPES-szel (Corning), kiegészítve 10% hő inaktivált magzati szarvasmarha szérummal (Sigma-Aldrich), 10 U/mL penicillin-streptomicinnel (Gibco) és 1× glutaMAX-szal (Gibco) párásított inkubátorban, 37 °C és 5% CO2 mellett. A sejteket rendszeresen teszteljük a mikoplazma jelenlétére a Lonza Walkersville MycoAlert Mycoplasma Detection Kit segítségével (utolsó teszt 2020 januárjában). A sejtvonalakat legfeljebb 20 passzázsig tartjuk kultúrában, és rövid tandemismétlődések profilozásával (JURKAT) vagy a TCRb-NOTCH1 transzlokáció kimutatására szolgáló PCR segítségével hitelesítjük (TCRBJ2S4CUTLL1F:5′-GGACCCCCGGCTCTCTCAGTGCT-3′, NOTCH1CUTTL1R:5′-TCCCGCCCTCCAAAATAAGG-3′). Az utolsó sejthitelesítésre 2020 februárjában került sor. A humán CD4+ T-sejteket az AllCells cégtől vásároltuk. A primer humán mintákat tájékozott beleegyezéssel gyűjtöttük, és a Padovai Egyetem Intézeti Felülvizsgálati Bizottsága, az Associazone Italiana di Ematologia e Oncologia Pediatrica és a Berlin-Frankfurt-Münster (AIEOP-BFM) ALL 2000/2006 gyermekgyógyászati klinikai vizsgálatok felügyelete mellett elemeztük. Minden betegtől a vizsgálatba való belépéskor a Helsinki Nyilatkozatnak megfelelően tájékozott beleegyezést kaptunk a megmaradt anyag kutatási célú felhasználásához.
Antitestek és reagensek
Western blotokat a következő antitestekkel végeztünk: Aktin és CTCF a Millipore Sigmától (C4 klón; 07-729) és hasított NOTCH1 (Val1744) a Cell Signaling Technologytól (4147). A ChIP-seq-et a következő antitestek felhasználásával végeztük: CTCF a Millipore Sigma-tól (07-729), H3K27Ac (8173S) és H3K27me3 (9733S) a Cell Signaling Technology-tól, valamint H3K4me1 (07-473) a Millipore-tól.
In situ Hi-C
Az in situ Hi-C-t CD4+ T sejteken, Jurkat, CUTLL1 és beteg xenograftokon végeztük a korábban leírtak szerint. Röviden, a sejteket 1%-os formaldehiddel térhálósítottuk 10 percig szobahőmérsékleten. Hi-C reakciónként 5 millió sejtet lizáltunk és a sejtmagokat permeabilizáltuk. A DNS-t a New England Biolabs-tól származó MboI-vel emésztettük (R0147M). Az emésztett fragmentumokat a Jena Bioscience (NU-835-BIO14-S) biotinilált d-ATP-vel jelöltük és ligáltuk. Az RNáz-kezelést és a keresztkötések visszafordítása érdekében végzett proteináz K kezelést követően a sejtmagokat Covaris E220 segítségével szonikáztuk, hogy a fragmentumok átlagos hossza 400 bp legyen. A Thermo Fisher Scientific (65001) Streptavidin gyöngyeit használtuk a biotinnal jelölt fragmentumok lehúzásához. A DNS tisztítását és izolálását követően a végleges könyvtárakat a NEBNext® Ultra™ II DNA Library Prep Kit for Illumina® segítségével készítettük el, és párosított végű szekvenálással szekvenáltuk 150 bp olvasási hosszúsággal egy Illumina HiSeq 2500-on, hogy mintánként átlagosan 400 millió olvasatot állítsunk elő.
ChIP-seq profilalkotás
CD4+ T-sejteket, Jurkat, CUTLL1 és beteg xenotranszplantátumokat 1% formaldehiddel és 1% magzati szarvasmarha-szérummal PBS-ben 10 percig szobahőmérsékleten keresztkötéssel kezeltük. A reakciót 0,2 M glicinnel oltottuk szobahőmérsékleten 5 percig. Ezután a sejteket PBS-szel mostuk és pelletáltuk.
A CTCF ChIP-ek esetében az immunprecipitációt a korábban leírt protokoll alapján végeztük. Az 50 millió sejtet tartalmazó pelletet 5 mL lízispufferrel (50 mM HEPES-KOH, pH 7,5, 140 mM NaCl, 1 mM EDTA, 10% glicerin, 0,5% NP-40, 0,25% Triton X-100) 10 percig 4 °C-on lizáltuk. A sejtmagokat 7 percig 1350×g-on pelletáltuk, majd 10 mM Tris pH 8, 1 mM EDTA és 0,1% SDS-ben reszuszpendáltuk. A kromatint a Covaris E220 rendszerrel 400 bp átlagos fragmentumhosszúra nyírtuk, majd 10 percig 15 000 rpm-en pörgettük az oldhatatlan kromatin és a törmelék eltávolítása érdekében. A felülúszót 20 μL Dynabeads Protein G-vel inkubáltuk 30 percig, mielőtt a gyöngyöket eldobtuk. A teljes térfogat egy százalékát mentettük el inputként, a többit pedig anti-CTCF antitesttel inkubáltuk egy éjszakán át. Összesen 100 μL Dynabeads Protein G-t adtunk 2 órára. A kötött fragmentumokat kétszer mostuk 1 mL alacsony sótartalmú pufferrel (20 mM Tris-HCl pH 8,0, 150 mM NaCl, 2 mM EDTA, 1% w/v Triton X-100, 0.1% w/v SDS), egyszer magas sótartalmú pufferrel (20 mM Tris-HCl pH 8,0, 500 mM NaCl, 2 mM EDTA, 1% w/v Triton X-100 és 0,1% w/v SDS), egyszer lítium-klorid pufferrel (10 mM Tris-HCl pH 8.0, 250 mM LiCl, 1 mM EDTA, 1% w/v NP-40 és 1% w/v dezoxikolsav), és kétszer TE-vel (10 mM Tris pH 8, 1 mM EDTA).
A hiszton ChIP-hez a sejteket 375 μL nukleáris inkubációs pufferben (15 mM Tris pH 7.5, 60 mM KCl, 150 mM NaCl, 15 mM MgCl2, 1 mM CaCl2, 250 mM szacharóz, 0,3% NP-40, 1 mM NaV, 1 mM NaF és 1 EDTA-mentes proteáz inhibitor tabletta (Roche)/10 mL H2O-ban) 10 percig jégen. A sejtmagokat egyszer mostuk emésztési pufferrel (10 mM NaCl, 10 mM Tris pH 7,5, 3 mM MgCl2, 1 mM CaCl2, 1 mM NaV, 1 mM NaF és 1 EDTA-mentes proteáz inhibitor tabletta (Roche)/10 mL H2O-ban), majd 57μL, 4,5 egység MNázt (USB) tartalmazó emésztési pufferben 1 órán át 37 °C-on reszuszpendáltuk. Az MNáz aktivitást 10 percig jégen oltottuk EDTA hozzáadásával 20 mM végső koncentrációig. A sejtmagokat pelletáltuk és 300μL Nuclei Lysis Bufferben (50 mM Tris-HCl pH 8,0, 10 mM EDTA pH 8,0, 1% SDS, 1 mM NaV, 1 mM NaF és 1 EDTA-mentes proteáz inhibitor tabletta (Roche)/10 mL H2O-ban) reszuszpendáltuk, majd 5 percig (30 s be, 30 s ki) Bioruptor Pico-val (Diagenode) szonikáltuk. A lizátumot 5 percig maximális sebességen centrifugáltuk a törmelék eltávolítása érdekében. A felülúszóhoz kilenc térfogat IP hígító puffert (0,01% SDS, 1,1% Triton X-100, 1,2 mM EDTA pH 8,0, 16,7 mM Tris-HCl pH 8,0, 167 mM NaCl, 1 mM NaV, 1 mM NaF és 1 EDTA-mentes proteáz inhibitor tabletta (Roche)/10 mL H2O-ban) adtunk. Összesen 50 μL Dynabeads Protein G-t adtunk hozzá, és a mintát 4 °C-on inkubáltuk 30 percig, forgatva. A minta egy százalékát megtartottuk bemenetként, a fennmaradó mintát pedig 3 csőre osztottuk. Összesen 50 μL Dynabeads Protein G-t adtunk 15 μL megfelelő ellenanyaggal konjugálva minden egyes csőbe, mielőtt egy éjszakán át 4 °C-on, forgatás mellett inkubáltuk volna. A gyöngyhöz kötött komplexeket egyenként 5 percig mostuk 1 mL alacsony sótartalmú pufferben, magas sótartalmú pufferben, LiCl pufferben és kétszer TE-vel.
A gyöngyhöz kötött komplexek eluálásához a gyöngyöket 50 μL elúciós pufferben (100 mM NaHCO3, 1% w/v SDS) reszuszpendáltuk, és 15 percig 65 °C-on inkubáltuk, 1000 RPM-es rázással termomixeren (Thermo Scientific). Az eluálást másodszor is megismételtük, majd 100 μL RNáz puffert (12 μL 5 M NaCl, 0,2 μL 30 mg/ml RNáz és 88 μL TE) adtunk minden egyes ChIP- és bemeneti mintához. A mintákat 37 °C-on inkubáltuk 20 percig, majd 100 μL proteináz K puffert (2,5 μL 20 mg/ml proteináz K, 5 μL 20% SDS és 92,5 μL TE) adtunk hozzá egy éjszakán át 65 °C-on. Ugyanekkora mennyiségű fenol:kloroform oldatot adtunk hozzá, és alaposan összekevertük. Az elegyet MaXtract High Density csövekbe (Qiagen) vittük át, és 8 percig 15 000 rpm-en centrifugáltuk. A felső fázist új csövekbe vittük át, és 1,5 μL 20 mg/mL glikogénnel, 30 μL 3M nátrium-acetáttal és 800 μL etanollal kevertük. A mintákat -80 °C-on inkubáltuk fagyasztásig, majd 4 °C-on 30 percig 15 000 rpm-en centrifugáltuk. A felülúszót eltávolítottuk, és a pelleteket 800 μL 70%-os jéghideg etanolban mostuk, majd 10 percig pörgettük 4 °C-on 15 000 rpm-en. Az etanol óvatos eltávolítása után a pelleteket levegőn szárítottuk, majd 30 μL 10 mM Tris pH 8-ban reszuszpendáltuk.
Az IP és az input DNS mennyiségét ezután Qubit 3.0 fluorométerrel számszerűsítettük. A könyvtárakat a KAPA HyperPrep Kit (KK8505) segítségével állítottuk elő és szekvenáltuk Illumina NextSeq 500 készülékkel, mintánként átlagosan 28 millió olvasás mélységig.
RNA-seq profilalkotás
Az RNS-t mintánként 3 millió sejtből izoláltuk a Bio-Rad Aurum™ Total RNA Mini Kit segítségével, és az Agilent RNA 6000 Nano Kit segítségével számszerűsítettük Agilent Bioanalyzerrel. A könyvtárakat rRNS-deplécióval állítottuk elő az Illumina TruSeq® Stranded mRNS Library Prep Kit segítségével alacsony koncentrációjú kiindulási minta esetén, és szekvenáltuk egyvégű szekvenálással Illumina NextSeq 500-on, mintánként átlagosan 18 millió olvasás mélységig.
DNS-metilációs profilalkotás
A genomi DNS-t az AllPrep DNA/RNS Micro Kit (Qiagen) segítségével izoláltuk. A genom-szintű DNS-metilációs állapot felméréséhez mRRBS-t végeztünk . A Qubit 3.0 műszerrel végzett fluorometriás mennyiségi meghatározást követően a genomi DNS-t az MspI restrikciós enzimmel (New England Biolabs) emésztettük, és a körülbelül 100-250 bázispár hosszúságú fragmentumokat szilárd fázisú reverzibilis immobilizációs (SPRI) gyöngyökkel (MagBio Genomics) méretszelektáltuk. Az így kapott DNS-t az EZ DNA Methylation-Lightning Kit (Zymo Research) segítségével biszulfitkonverziónak vetettük alá. A biszulfit-konvertált egyszálú DNS-ből a Pico Methyl-Seq Library Prep Kit (Zymo Research) segítségével könyvtárakat hoztunk létre, amelyeket ezután összevontunk az Illumina NextSeq 500 műszeren történő szekvenáláshoz a NextSeq 500/550 V2 High Output reagens kit használatával (1 × 75 ciklus), mintánként legalább 50 millió leolvasási mélységig.
Teljes genom szekvenálás
A sejtvonalakból vagy betegmintákból származó hárommillió sejtet pelletáltuk, és 1 ml Cell Lysis Solution (Qiagen) oldatban (500 μg RNáz A-val keverve) reszuszpendáltuk. 37 °C-on 15 percig végeztük a lízisreakciót. Összesen 333 μL fehérjecsapadék-leválasztó oldatot (Qiagen) adtunk minden egyes mintához, amelyet vortexeltünk, majd 10 percig 2000×g-en centrifugáltunk. A felülúszót 1 mL izopropanollal kevertük, amíg a DNS-szálak ki nem csapódtak az oldatból. A felülúszó elvetése után a DNS-pelletet 1 mL 70%-os etanollal mostuk, majd 1 percig 2000×g-nél centrifugáltuk. Ezután az etanolt kiöntöttük, és a pelletet 15 percig levegőn szárítottuk, majd 50-100 μL DNS-hidratáló oldatban (Qiagen) reszuszpendáltuk. A DNS-t párosított Illumina szekvenálással szekvenáltuk 30×-os lefedettséggel.
Immunoprecipitáció
Minden immunprecipitációs reakcióhoz összesen 100 millió sejtet pelletáltunk és 10 percig A pufferben (10 mM HEPES pH 8,0, 1,5 mM MgCl2, 10 mM KCl, 0,5 mM DTT) inkubáltuk jégen. Ezután a sejteket 12 ütés hatására lizáltuk egy 7 ml-es laza pisztolyos szövetdarálóval (Wheaton, 357542) és 7 percig 2000 rpm-en centrifugáltuk. A sejtmag pelleteket 5 térfogat TENT pufferben (50 mM Tris pH 7,5, 5 mM EDTA, 150 mM NaCl, 1% Triton X-100, 5 mM MgCl2) reszuszpendáltuk, és benzonázzal kezeltük 30 percig, mielőtt 5 átjárást végeztünk 25 g × 5/8 in. fecskendőn keresztül. Az oldhatatlan frakciót 2000 rpm-en 7 percig tartó centrifugálást követően eltávolítottuk, és egy éjszakán át antitesttel hibridizált Dynabeads Protein G-vel inkubáltuk. Összesen 2 millió sejtet távolítottunk el a bemenethez. A gyöngyöket és a sejtmag-lizátumokat 6-szor mostuk TENT pufferrel, majd 0,1 M glicin pH 2,5-ben eluáltuk 100 mM Tris pH 8,0-val előtte. NuPAGE LDS mintapuffert adtunk az eluátumokhoz és a bemenetekhez, amelyeket aztán 15 percig inkubáltunk 70 °C-on, mielőtt western blottal elemeztük őket.
Nyilvános adatgyűjtés
A nyilvános CTCF ChIP-seq adatokat a Cistrome Data Browserből (a csúcsfájlokhoz) és az NCBI GEO-ból (a fastq-fájlokhoz, Additional file 2: S1 táblázat) gyűjtöttük. A hiszton-modifikációs ChIP-seq adatokat az NCBI GEO-ból és az ENCODE-ból gyűjtöttük (a bam-fájlokhoz). A több sejttípusra vonatkozó nyilvános RNS-seq adatokat az ENCODE-ból gyűjtöttük (fastq fájlok esetében). A DNS-metilációs profilalkotási adatokat az ENCODE-ból (a bed bedMethyl-fájlokhoz) és az NCBI GEO-ból gyűjtöttük. A Hi-C adatokat az NCBI GEO-ból és az ENCODE-ból gyűjtöttük (fastq fájlok esetében). Az ATAC-seq adatokat az NCBI GEO-ból gyűjtöttük (fastq fájlok esetében). A BRCA, COAD, LUAD és PRAD minták teljes genomszekvenálási adatait az International Cancer Genome Consortium (ICGC) adatportáljáról gyűjtöttük. Részletes információk, beleértve az e munka során gyűjtött összes nyilvános adatkészlet hozzáférési azonosítóját, a 6. kiegészítő fájlban találhatók: S5. táblázat.
Adatok feldolgozása
ChIP-seq adatok elemzése
A fastq fájlokban lévő ChIP-seq adatok szekvenciaillesztését a Cistrome DB-ben használt standard elemzési csővezetékkel végeztük el a konzisztencia és a reprodukálhatóság érdekében. Minden szekvenciaadatok genomiális összehangolását a Chilin pipeline segítségével végeztük el, alapértelmezett paraméterekkel ($ chilin simple -p narrow -s hg38 –threads 8 -t IN.fq -i PRENAME -o OUTDIR). Röviden, a szekvenciaolvasásokat a humán referencia genomhoz (GRCH38/hg38) igazítottuk a BWA segítségével ($ bwa aln -q 5 -l 32 -k 2 -t 8 INDEX IN.fq > PRENAME.sai $ bwa {samse | sampe} INDEX PRENAME.sai IN.fq > PRENAME.sam). A Sam-fájlokat ezután a samtools segítségével bam-fájlokká alakítottuk át ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). A CTCF ChIP-seq adatkészletek esetében a MACS2-t használtuk a csúcsok hívására a 0,01-es FDR küszöbérték alatt ($ macs2 callpeak –SPMR -B -q 0,01 –keep-dup 1 -g hs -t PRENAME.bam -n PRENAME –outidr OUTDIR). A legalább 4 szeres feldúsulással rendelkező csúcsokat megtartottuk. A Bigwiggle fájlokat a BEDTools és az UCSC eszközeivel generáltuk ($ bedtools slop -i PRENAME.bdg -g CHROMSIZE -b 0|bedClip stdin CHROMSIZE PRENAME.bdg.clip $ LC_COLLATE=C sort -k1,1 -k2,2n PRENAME.bdg.clip > PRENAME.bdg.sort.clip $ bedGraphToBigWig PRENAME.bdg.sort.clip CHROMSIZE PRENAME.bw). Végül csak a legalább 2000 csúcsot tartalmazó CTCF ChIP-seq minták kerültek be a downstream integratív elemzésbe.
ATAC-seq adatelemzés
A nyers szekvenálási olvasatok trimmelésére a Trim Galore-t használtuk ($ trim_galore –nextera –phred33 –fastqc –paired R1.fq R2.fq -o OUTDIR). A leolvasásokat a Bowtie2 segítségével igazítottuk a humán referencia genomhoz (GRCH38/hg38) ($ bowtie2 -p 10 -X 2000 -x INDEX -1 R1.fq -2 R2.fq -S PRENAME.sam). A Sam-fájlokat ezután samtools segítségével bam-fájlokká alakítottuk át ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). Bedtools segítségével a bam fájlokat ágy formátumba konvertáltuk ($ bamToBed -i PRENAME.bam -bedpe > PRENAME_PE.bed). A mitokondriumok DNS-éhez leképezett leolvasásokat a downstream elemzésből kizártuk.
RNS-seq adatok elemzése
Az RNS-seq adatkészleteket a Salmon segítségével dolgoztuk fel ($ salmon quant –gcBias -i INDEX -l A -p 8 {-1 R1.fq -2 R2.fq| -r IN.fq} -o OUTDIR). A transzkriptom indexet a humán referencia genomra (GRCH38/hg38) építettük fel. A transzkript-szintű gyakorisági becsléseket génszintre összegeztük a differenciális expresszióelemzéshez használt “tximport” csomag segítségével. A DESeq2-t használtuk a differenciálisan expresszált gének azonosítására, és a különböző elemzésekben használt különböző küszöbértékeket ennek megfelelően tüntettük fel a kéziratban.
Hi-C adatelemzés
A Hi-C adatokat a HiC-Pro programmal dolgoztuk fel ($ HiC-Pro -i INDIR -o OUTDIR -c CONFIG -p). A kontakttérképeket 5 kb felbontással hoztuk létre. A nyers mátrixadatokat a Kromatin kölcsönhatások normalizálása című fejezetben leírt megközelítéssel normalizáltuk.
DNS-metilációs adatok elemzése
A DNS-metilációs adatokat (T-ALL sejtvonalak és T-ALL betegek esetében) bcl2fastq segítségével demultiplexáltuk, majd az 5′ végéről 10 bázispárt trimmeltünk a primer és adaptor szekvenciák eltávolítása érdekében TrimGalore segítségével. A szekvenciaillesztést a GRCh38/hg38 referencia genomhoz és a metilációs hívásokat a Bismark programmal végeztük ($ bismark –multicore 8 –bowtie2 -q -N 1 INDEX INFILE.fq). A CpG-kontextusban lévő citozinok lefedettségi (counts) fájljait Bismark segítségével generáltuk ($ bismark_methylation_extractor –multicore 8 –comprehensive –bedGraph INFILE_bismark_bt2.bam).
A teljes genom szekvenálási adatok elemzése
A teljes genom szekvenálási adatokból két T-ALL sejtvonal (Jurkat és CUTLL1) és két T-ALL beteg minta esetében azonosítottunk mutációkat. Az Illumina rövid leolvasású szekvenciáit a humán referencia genomhoz (GRCH38/hg38) igazítottuk BWA mem segítségével. SAMBlaster segítségével azonosítottuk a diszkordáns párokat, felosztottuk a leolvasásokat, és megjelöltük a feltételezett PCR-duplikátumokat. A SAMBAMBA-t használtuk a SAM igazítottak BAM formátumba való konvertálásához, és a samtools segítségével rendeztük az igazítottakat, hogy létrehozzuk az egyes mintáknak megfelelő BAM fájlt.
A VarDict segítségével azonosítottuk azokat a variánsokat, amelyek átfedésben voltak az unió CTCF-kötőhelyeivel. Az összes alapértelmezett paramétert használtuk, kivéve a “-f 0.1” paramétert, amelyet az olyan variánsok azonosítására használtunk, amelyeket az adott helyen az olvasások több mint 10%-a támogatott. A variánsokat a Variant Effect Predictor (VEP) segítségével annotáltuk, és egyéni szkripteket használtunk a TF-kötődést befolyásoló variánsok azonosítására.
Ismét a VarDict-et használtuk a CTCF- és NOTCH1-gének variánsainak azonosítására a négy minta esetében. Az összes alapértelmezett paramétert használtuk, kivéve a “-f 0.1” paramétert, amelyet az adott helyen a leolvasások több mint 10%-a által támogatott változatok azonosítására használtunk. A variánsokat a Variant Effect Predictor (VEP) segítségével annotáltuk, majd megszűrtük, hogy azonosítsuk azokat a mutációkat, amelyek (a) nem fordultak elő több mint 1%-ban egyetlen normális emberi populációban sem, vagy (b) a deleteriabilitás CADD pontszáma > 20 volt, vagy (c) jelen voltak a COSMIC adatbázisban.
Integratív modellezés és statisztikai elemzés
A CTCF-kötési repertoár azonosítása a humán genomban
A CTCF ChIP-seq számára összesen 793 adatkészletet gyűjtöttünk össze, köztük 787 nyilvános és 6 általunk generált adatkészletet (Additional file 2: S1 táblázat). Összesen 771, 2000-nél több csúcsot tartalmazó CTCF ChIP-seq adatkészletet használtunk fel ebben a vizsgálatban. Minden egyes adatkészlet 2050 és 198 021 közötti tartományban adhat MACS2-azonosított CTCF-csúcsokat, a medián 46 451 és összesen 36 873 077 csúcsot (Additional file 1: S1a ábra). A 771 adatkészletből származó 36 873 077 csúcs szomszédos CTCF-csúcsok csúcsai közötti intervallumhosszúságok eloszlása egy inflexiós pontot mutat ~ 150 bp-nél (Additional file 1: S1c ábra), ami az azonos kötőhely és a különböző kötőhelyek közötti határt jelzi. Ezért a CTCF-csúcsok összevonásához 150 bps-t használtunk határértékként. A gyakorlatban az egyes csúcsok csúcsától ± 75 bps-t kiterjesztettünk, hogy létrehozzunk egy 150 bp-os régiót, amelynek középpontja a csúcson van, és minden átfedő csúcs régiót egyesítettünk, hogy létrehozzuk a CTCF kötőhelyek egyesített halmazát, amely 688 429 nem átfedő helyet tartalmaz. Minden kötőhelyhez egy CTCF-foglaltsági pontszámot rendeltünk, amelyet azon ChIP-seq-adatkészletek számaként határoztunk meg, amelyek az adott helyen belül csúcsot mutatnak. Ennek megfelelően a foglaltsági gyakoriságot a foglaltsági pontszám és a CTCF ChIP-seq adatkészletek teljes számának hányadosaként határoztuk meg. Az azonosított CTCF-kötőhelyek robusztusságának további biztosítása érdekében 285 467 ≥ 3 foglaltsági pontszámú, nagy megbízhatóságú helyet választottunk ki a downstream elemzésekhez. A szövetségi kötőhelyeken belüli CTCF-motívumokat FIMO-val kerestük Jaspar mátrixszal (ID: MA0139.1), 1e-4 p-érték küszöbértékkel. Minden egyes CTCF-kötőhelyhez a legkisebb p-értékkel rendelkező motívumot megtartottuk.
Konstitutív CTCF-kötőhelyek azonosítása
A 285 467 CTCF-kötőhely foglaltsági pontszámának eloszlása (Additional file 1: S1d ábra, kék görbe) azt mutatja, hogy a CTCF-kötőhelyek többsége csak néhány adatsorban fordul elő, és a foglaltsági pontszám növekedésével csökken a kötőhelyek száma, ha a foglaltsági pontszám kicsi. Vannak azonban olyan CTCF-kötőhelyek, amelyek szinte minden adatkészletben erősen konzerváltak (pl. 600-nál nagyobb foglaltsági pontszámú kötőhelyek). Az 1. kiegészítő fájlban látható eloszlási görbe (kék) illesztéséhez hatványtörvény-függvényt használunk: A konstitutív CTCF-helyek határértékének meghatározásához. Az Oi-t az i-nek megfelelő foglaltsági pontszámú megfigyelt CTCF-kötőhelyek számaként jelöljük, az Ei-t pedig az i-nek megfelelő foglaltsági pontszámú várható CTCF-helyek számaként. Az Oi adatokra való hatványtörvény-illesztés a következőképpen írható le (Additional file 1: S1d ábra, zöld):
A konstitutív CTCF-kötőhelyek A határértékét a következőképpen határozzuk meg:
Más szóval, az A-nál nagyobb foglaltsági pontszámmal rendelkező összes megfigyelt CTCF-helynek hatszor többnek kell lennie a vártnál. Ezután A = 615 értéket határoztunk meg, és 80%-os foglaltsági gyakorisági határértéket használtunk a 22 097 konstitutív CTCF-kötőhely meghatározásához, ami megfelel az ≥ 616-os foglaltsági pontszámnak mind a 771 CTCF ChIP-seq adatkészletben.
Rák-specifikus nyert/vesztett CTCF-kötőhelyek azonosítása
A következő 2 kritériumot használtuk a rák-specifikus vesztett CTCF-kötőhelyek azonosításához: (1) A CTCF-kötőhelynek alacsonyabb foglaltsági gyakorisággal kell rendelkeznie az adott ráktípushoz tartozó adatkészletekben az összes adatkészlet foglaltsági gyakoriságához képest, és (2) a CTCF-kötési szint (normalizált ChIP-seq olvasásszámként számszerűsítve) az adott helyen alacsonyabb a rákos adatkészletekben, mint a többi adatkészletben. A nyert CTCF-helyek esetében a kritériumok fordítottját használtuk. Röviden, minden rákos típusban minden egyes CTCF-kötőhelyre kiszámítottuk a rákos adatkészletekben elért foglaltsági pontszámot, valamint az összes 771 adatkészletben elért foglaltsági pontszámot. A CTCF-kötési szinteket egy normalizált olvasásszám-mátrixból kaptuk, amelyben először a ChIP-seq olvasásszámokat (RPKM) számoltuk ki az összes adatkészletben lévő unió CTCF-kötőhelyekre, majd ezt követően kvantilis normalizálást végeztünk. Párosítatlan kétfarkú Student’s t-tesztet használtunk a kötődési szintek különbségének számszerűsítésére az adatkészletek különböző csoportjai között, majd a p-értéket a Benjamini-Hochberg-eljárás segítségével korrigáltuk . Ezenkívül a kötési foglaltsági pontszámokat és a kötési szinteket összehasonlítottuk a rákos és a megfelelő normál szövetekből vagy sejttípusokból származó adatkészletek között, hogy figyelembe vegyük a szövetspecifikusság, nem pedig a rákspecifikusság potenciális zavaró tényezőjét. A rák-specifikus CTCF-kötőhelyek azonosításának részletes kritériumait az alábbiakban ismertetjük:
-
Rák-specifikus elveszett CTCF-kötőhelyek: (1) foglaltsági gyakoriság ≤ 0,2 a rákos adatkészletekben; (2) foglaltsági gyakoriság ≥ 0,7 a 771 adatkészletben; (3) foglaltsági gyakoriság ≥ 0.5 (foglaltsági pontszám ≥ 2) a megfelelő normál szöveti adatkészletekben; (4) a CTCF-szintek alacsonyabbak a rákban az összes többi adatkészlethez képest (statisztikai pontszám < 0), (5) a CTCF-szintek alacsonyabbak a rákban a megfelelő normál szöveti adatkészletekhez képest (statisztikai pontszám < 0), (6) az átlagos CTCF-kötési jelek (RPKM) < 5 a rákos adatkészletekben.
-
Rák-specifikusan szerzett CTCF-kötőhelyek: (1) foglaltsági gyakoriság ≥ 0,5 (foglaltsági pontszámmal ≥ 2) a rákos adatkészletekben, (2) foglaltsági gyakoriság ≤ 0,2 a 771 adatkészletben, (3) foglaltsági pontszám = 0 a megfelelő normál szöveti adatkészletekben, (4) a CTCF-szintek szignifikánsan magasabbak a rákban az összes többi adatkészlethez képest (FDR ≤ 0.01), (5) a CTCF kötődési szintek szignifikánsan magasabbak a rákban az illesztett normál szöveti adatkészletekhez képest (FDR ≤ 0,01), (6) az átlagolt CTCF kötődési jelek (RPKM) > 2 a rákos adatkészletekben.
A specifikusan nyert és elveszett CTCF-kötőhelyeket az egyes ráktípusok esetében a 4. kiegészítő fájl: S3. táblázat tartalmazza.
A differenciális kromatin-hozzáférhetőség kvantitatív meghatározása
A CTCF-kötőhelyek körüli differenciális kromatin-hozzáférhetőség értékeléséhez a Ref. feldolgozott adatait használtuk, amelyek a TCGA pán-rákcsúcskészleten belüli normalizált ATAC-seq inszerciós számok mátrixát tartalmazzák. A BRCA, CRC, LUAD és PRAD egyes ráktípusai közül a pánrák ATAC-seq-csúcsokat, amelyek átfedésben vannak az azonosított ráktípus-specifikus elveszett vagy nyert CTCF-kötőhelyekkel, használtuk a downstream elemzésekhez. Az egyes csúcsok ATAC-seq differenciális pontszámát a megfelelő ráktípusba tartozó betegmintákból származó, normalizált ATAC-seq inszerciós számok átlagának a más ráktípusokba tartozó betegmintákhoz képest bekövetkezett hajtott változásaként számszerűsítettük, majd az ATAC-seq differenciális pontszámot a CTCF-kötőhelyekkel átfedő csúcshoz rendeltük.
A következetesség érdekében a TCGA ATAC-seq-adatokhoz használt megközelítést alkalmaztuk a T-ALL sejtvonal Jurkat és normál CD4+ T-sejtekből gyűjtött ATAC-seq-adatok elemzésére is. A Jurkat és a T-sejtek összes adatkészletére vonatkozóan adatmátrixot hoztunk létre az ATAC-seq nyers olvasatszámának felhasználásával az unió CTCF-kötőhelyein. A log2 skálázott mátrixon kvantilis normalizálást alkalmaztunk (álszám = 5). Az ATAC-seq differenciális pontszámot a Jurkat versus CD4+ T-sejt adatkészletek közötti átlagolt normalizált ATAC-seq-számok szeres változásaként mértük az egyes CTCF-kötőhelyeken.
Kromatin kölcsönhatások normalizálása
Adott egy Hi-C kontakt térkép A = {aij}, az aij pontszám a két genomi régió i és j közötti leképezett olvasásokat tükrözi. Tegyük fel, hogy a bin méret 5 kb, az i és j régiók genomi távolsága ∣i – j ∣ × 5kb. Mivel a két bin közötti érintkezési valószínűség a genomiális távolság növekedésével csökken, a kontakttérképet a következőképpen normalizáltuk: bármely adott dk = k × 5kb genomiális távolság esetén egy \( {\overline{S}}_{d_k} \) normalizációs tényezőt számszerűsítünk, mint az ugyanazon kromoszómában lévő, azonos dk genomiális távolsággal rendelkező összes bin-pár közötti kölcsönhatások átlagát, például, \( {\overline{S}}}_{d_k}=\left({\sum}_{j-i=k}{a}_{ij}\right)/n \), ahol n a dk távolságú bin-párok teljes száma. A két dk távolságú bins közötti aij kölcsönhatási pontszámot ezután \( {\overline{S}}}_{d_k} \) normáltuk a következő módon: \( {a}_{ij}^{\prime }={a}_{ij}/{\overline{S}}}_{d_k} \). Ezt a megközelítést alkalmazva az A mátrixot minden kromoszómán belül \( A^{\prime }=\left\{{a}_{ij}^{\prime}\right\} \) mátrixra normalizáltuk.
A differenciális kromatin kölcsönhatások kimutatása
A rákos és a normál adathalmazban a normalizált Hi-C kontakt térképeket C = {cij} és N = {nij} jelekkel jelöltük. Egy adott x CTCF-kötőhelyhez (xc koordinátával) és egy előre meghatározott L genomiális távolsághoz az x és a közeli, L-ig terjedő genomiális távolsággal rendelkező, nem átfedő 5 kb-os tartományok közötti kromatin kölcsönhatásokat a C és az N tartományból gyűjtjük össze. Pontosabban, az x és a C-ben lévő közeli 5 kb-os bins közötti kölcsönhatási pontszámokat a következőképpen gyűjtjük össze: IC = {cij} , míg i vagy j egyenlő ⌊xc/5kb⌋, és 0 < (j – i) × 5kb ≤ L. Hasonlóképpen, az x és az N-ben lévő közeli 5-kb-os bins közötti kölcsönhatási pontszámokat IN = {nij} formában gyűjtöttük. Ezután párosított kétfarkú Student’s t-tesztet alkalmaztunk az IC-re és az IN-re, hogy számszerűsítsük az x CTCF-kötőhelyet körülvevő rákos és normál sejtek közötti eltérő kölcsönhatást.
A CTCF-kötés és a génexpresszió kapcsolata
A CTCF-kötés és a génexpresszió közötti kapcsolat vizsgálatához összesen 54 olyan sejttípust választottunk ki (Additional file 6: Table S5), amelyek esetében mind a CTCF ChIP-seq-adatok, mind az RNS-seq-adatok nyilvánosan elérhetőek (Additional file 6: Table S5) az egyes CTCF-génpárok esetében ugyanazon kromoszómán belül. A CTCF-kötődés szintjének meghatározásához a ChIP-seq-adatokból származó egyesülési CTCF-kötőhelyeken a kilobázismilliónkénti olvasásszám-mátrixot (RPKM) használtuk. A leolvasásszám-mátrixot az RPKM négyzetgyökével skáláztuk, majd kvantilis normalizálást végeztünk. A génexpressziós szintet minden gén esetében az RNA-seq-adatokból származó transzkriptumok millióra jutó négyzetgyökének (TPM) felhasználásával mértük. Minden CTCF-gén pár esetében a CTCF-hely és a gén közötti asszociációt mind az 54 sejttípusra vonatkozóan a normalizált CTCF-kötési szint és a génexpresszió közötti R korrelációs együttható segítségével számszerűsítettük (3a. ábra). A CTCF-gén párokat “erősen korreláltnak” tekintettük, ha az R2 nagyobb volt, mint 0,25, pl, korrelációs együttható nagyobb, mint 0,5 vagy kisebb, mint – 0,5, és a magasan korrelált CTCF-gén párok az összes CTCF-gén pár 1,3%-át teszik ki (Additional file 1: Fig. S8a).
Konstitutív CTCF-kötőhelyekkel határolt kromatindomének azonosítása
Minden CTCF-kötőhely esetében a hozzá tartozó kromatindomént úgy definiáltuk, mint azt a genomi régiót, amely (1) tartalmazza az adott CTCF-kötőhelyet, (2) ellentétes irányú motívumokkal rendelkező konstitutív CTCF-kötőhelyek párja határolja, és (3) legalább 100 kb és legfeljebb 1 MB területet foglal el a CTCF-kötőhely mindkét oldalán. A 3b. ábra tartalmazza a konstitutív CTCF-kötőhelyekkel határolt kromatindomének meghatározásának sematikus ábráját.
A CTCF-kötőhelyeket körülvevő DNS-metilációs változások kimutatása
A DNS-metilációs változásokat az egyes CTCF-kötőhelyekre központosított 300 bp hosszúságú régióban detektáltuk. Mind a rákos sejtvonalakban, mind a megfelelő normál szövetekben megtartottuk azokat a régiókat, amelyekben legalább 3 CpG-t legalább 5 leolvasás (≥ 5×) fedett le. Egy 300 bp hosszúságú régiót akkor detektáltunk differenciálisan metiláltnak, ha a régióban található összes CpG (≥ 5×) átlagos differenciális metilációs szintje nagyobb volt, mint 20% .
Mutációs arány és differenciális motívum pontszám meghatározása
Minden CTCF-kötőhely esetében a nyers mutációs számot az összes mintában/betegben a mutációs események előfordulása alapján számoltuk ki a CTCF-kötőhelyre központosított 400 bp hosszúságú régióban található minden egyes bázispárban. A CTCF-kötőhelyek egy csoportjára vonatkozó mutációs rátát a CTCF-kötőhelyek számának átlagolt mutációs számaként számoltuk ki a 400 bp-os régión belüli minden egyes bázispárra vonatkozóan.
A motívum pontszámát a CTCF-pozíció súlymátrixának (Jaspar , Matrix ID: MA0139.1) a CTCF-motívum vagy a CTCF-kötőhely középpontjában lévő 19 bp-os DNS-szekvenciára történő pontozásával, log valószínűségi arányok segítségével (a háttér nukleotidok gyakorisága az A,C,G,T-hez hasonlóan). A differenciális motívumpontszámot a referencia és a mutált szekvenciák motívumpontszámainak összehasonlításával számoltuk ki.
DNS-szekvencia motívumelemzés
A DNS-szekvencia motívumgazdagodási elemzését az MDSeqPos (1.0.0 verzió) segítségével végeztük a Cistrome-on, alapértelmezett paraméterekkel (-cisrome -Homo Sapien vagy Mus musculus). A de novo motívumelemzéseket a HOMER (4.10-es verzió) findmotifs.pl moduljával és a MEME (5.1. verzió.1) a következő paraméterekkel: meme -dna -mod zoops -maxw 20 -evt -0.01.
A CTCF domainon belüli differenciálisan kölcsönhatásban lévő régiók azonosítása
A CTCF-kötőhelyek adott készletére vonatkozóan összegyűjtöttük a CTCF-kötőhelyek és az egyes domainon belüli, nem átfedő bins közötti kromatin kölcsönhatási változásokat, amelyeket a rákos sejtek normalizált Hi-C kontakt-térképéből mértünk a megfelelő normál sejtekhez képest (Additional file 1: Fig. S14b). A rák-specifikus elveszett CTCF-kötőhelyekkel csökkent kölcsönhatásokat (log2 FC < -1, átlagos log2 kölcsönhatás > 0) mutató régiókat és a rák-specifikus nyert CTCF-kötőhelyekkel fokozott kölcsönhatásokat (log2 FC > 1, átlagos log2 kölcsönhatás > 0) mutató régiókat használtuk a downstream transzkripciós faktor (TF) dúsítási elemzéshez.
Transzkripciós faktorok dúsítási elemzése
A TF-dúsítási elemzéshez a BART algoritmus felülvizsgált változatát használtuk. Röviden, az unió DNase I hiperszenzitív helyeinek (UDHS) gyűjteményét korábban a humán genom összes jelölt cisz-szabályozó elemének repertoárjaként kuratálták, és 7032 ChIP-seq adatkészletet gyűjtöttek 883 TF-hez , minden egyes TF-hez egy vagy több ChIP-seq adatkészletet gyűjtöttek több sejttípusból vagy állapotból. Minden TF számára bináris profilt hoztunk létre az UDHS-en, amely jelzi, hogy a TF bármelyik ChIP-seq adatkészletéből legalább egy csúcs található-e az egyes UDHS-en belül. A kötődési gazdagodási elemzést minden TF esetében úgy alkalmaztuk, hogy összehasonlítottuk a TF kötődését a kiválasztott genomi régiókkal átfedő UDHS-ek egy részhalmazán a TF kötődésével az UDHS-en. p-értéket kaptunk a kétfarkú Fisher-féle egzakt teszt segítségével.