Experimentele procedure
Patient xenografting and cell culture
De menselijke T-ALL cellijnen omvatten CUTLL1 (gift van Adolfo Ferrando, Columbia University) en JURKAT (American Type Culture Collection (ATCC), Manassas, VA, #CCL-119) . De cellen werden gekweekt in RPMI1640-medium met l-glutamine en 25 mM HEPES (Corning), aangevuld met 10% door warmte geïnactiveerd foetaal runderserum (Sigma-Aldrich), 10 U/mL penicilline-streptomycine (Gibco), en 1× glutaMAX (Gibco) in een bevochtigde incubator bij 37 °C en 5% CO2. De cellen worden periodiek getest op de aanwezigheid van mycoplasma met de Lonza Walkersville MycoAlert Mycoplasma Detection Kit (laatste test in januari 2020). De cellijnen worden maximaal 20 passages in cultuur gehouden en worden geauthenticeerd met behulp van short-tandem repeats profiling (JURKAT) of met behulp van PCR om de TCRb-NOTCH1 translocatie (TCRBJ2S4CUTLL1F:5′-GGACCCGGCTCTCAGTGCT-3′, NOTCH1CUTTL1R:5′-TCCCGCCCTCCAAAATAAGG-3′) op te sporen. De laatste celverificatie werd uitgevoerd in februari 2020. Menselijke CD4+ T-cellen werden gekocht bij AllCells. Primaire menselijke monsters werden verzameld met geïnformeerde toestemming en geanalyseerd onder toezicht van de Institutional Review Board van de Universiteit van Padua, de Associazone Italiana di Ematologia e Oncologia Pediatrica, en de Berlijn-Frankfurt-Münster (AIEOP-BFM) ALL 2000/2006 pediatrische klinische proeven. Geïnformeerde toestemming voor het gebruik van restmateriaal voor onderzoeksdoeleinden werd verkregen van alle patiënten bij het begin van de proef, in overeenstemming met de Verklaring van Helsinki.
Antilichamen en reagentia
Western blots werden uitgevoerd met de volgende antilichamen: Actine en CTCF van Millipore Sigma (kloon C4; 07-729) en gekloofde NOTCH1 (Val1744) van Cell Signaling Technology (4147). ChIP-seq werden uitgevoerd met de volgende antilichamen: CTCF van Millipore Sigma (07-729), H3K27Ac (8173S), en H3K27me3 (9733S) van Cell Signaling Technology, en H3K4me1 (07-473) van Millipore.
In situ Hi-C
In situ Hi-C werd uitgevoerd op CD4 + T-cellen, Jurkat, CUTLL1, en patiënt xenograften zoals eerder beschreven . In het kort, werden cellen verknoopt met 1% formaldehyde gedurende 10 min bij kamertemperatuur. Per Hi-C reactie werden 5 miljoen cellen gelyseerd en de kernen permeabiliseerd. Het DNA werd gedigesteerd met MboI van New England Biolabs (R0147M). De verteerde fragmenten werden gelabeld met gebiotinyleerd d-ATP van Jena Bioscience (NU-835-BIO14-S) en geligeerd. Na RNase behandeling en Proteinase K behandeling om crosslinks om te keren, werden kernen gesonificeerd met behulp van een Covaris E220 om een gemiddelde fragmentlengte van 400 bp te produceren. Streptavidine kralen van Thermo Fisher Scientific (65001) werden gebruikt om naar beneden te trekken biotine gelabelde fragmenten. Na zuivering en isolatie van DNA werden de uiteindelijke bibliotheken bereid met de NEBNext® Ultra™ II DNA Library Prep Kit voor Illumina® en gesequeneerd via gepaarde end sequencing met een leeslengte van 150 bp op een Illumina HiSeq 2500 om gemiddeld 400 miljoen lezingen per monster te produceren.
ChIP-seq profilering
CD4 + T-cellen, Jurkat, CUTLL1, en patiënt xenograften werden crosslinked met 1% formaldehyde en 1% foetaal runderserum in PBS gedurende 10 min bij kamertemperatuur. De reactie werd gedoofd met 0,2 M glycine bij kamertemperatuur gedurende 5 min. Cellen werden vervolgens gewassen met PBS en gepelletiseerd.
Voor CTCF ChIPs, immunoprecipitatie werd uitgevoerd op basis van een protocol eerder beschreven. Een pellet met 50 miljoen cellen werd gelyseerd met 5 mL lysisbuffer (50 mM HEPES-KOH, pH 7,5, 140 mM NaCl, 1 mM EDTA, 10% glycerol, 0,5% NP-40, 0,25% Triton X-100) gedurende 10 min bij 4 ° C. Nuclei werden gepelletiseerd bij 1350×g gedurende 7 min en geresuspendeerd in 10 mM Tris pH 8, 1 mM EDTA, en 0,1% SDS. Chromatine werd geschoren met een Covaris E220 systeem tot een gemiddelde fragmentlengte van 400 bp en gesponnen bij 15.000 rpm gedurende 10 min om onoplosbare chromatine en puin te verwijderen. Het supernatant werd geïncubeerd met 20 pi van Dynabeads Protein G gedurende 30 minuten voor het weggooien van de kralen. Een procent van het totale volume werd opgeslagen als input en de rest werd geïncubeerd met anti-CTCF antilichaam overnacht. In totaal werd 100 ul Dynabeads Protein G toegevoegd gedurende 2 uur. Gebonden fragmenten werden tweemaal gewassen met 1 ml zoutarme buffer (20 mM Tris-HCl pH 8,0, 150 mM NaCl, 2 mM EDTA, 1% w/v Triton X-100, en 0.1% w/v SDS), eenmaal met hoogzoutbuffer (20 mM Tris-HCl pH 8.0, 500 mM NaCl, 2 mM EDTA, 1% w/v Triton X-100, en 0,1% w/v SDS), eenmaal met lithiumchloridebuffer (10 mM Tris-HCl pH 8.0, 250 mM LiCl, 0,1% w/v Triton X-100, en 0,1% w/v SDS).0, 250 mM LiCl, 1 mM EDTA, 1% w/v NP-40, en 1% w/v deoxycholzuur), en tweemaal met TE (10 mM Tris pH 8, 1 mM EDTA).
Voor histon ChIPs werden de cellen gelyseerd in 375 μL nuclei incubatiebuffer (15 mM Tris pH 7.5, 60 mM KCl, 150 mM NaCl, 15 mM MgCl2, 1 mM CaCl2, 250 mM sacharose, 0,3% NP-40, 1 mM NaV, 1 mM NaF, en 1 EDTA-vrije protease inhibitor tablet (Roche)/10 mL in H2O) gedurende 10 min op ijs. Nuclei werden eenmaal gewassen met digest buffer (10 mM NaCl, 10 mM Tris pH 7,5, 3 mM MgCl2, 1 mM CaCl2, 1 mM NaV, 1 mM NaF, en 1 EDTA-vrije protease inhibitor tablet (Roche) / 10 mL in H 2 O) en geresuspendeerd in 57-uryl Digest Buffer met 4,5 eenheden MNase (USB) gedurende 1 uur bij 37 ° C. MNase-activiteit werd gedoofd gedurende 10 min op ijs na de toevoeging van EDTA tot een eindconcentratie van 20 mM. Nuclei werden gepelleteerd en geresuspendeerd in 300-uril Nuclei Lysis Buffer (50 mM Tris-HCl pH 8,0, 10 mM EDTA pH 8,0, 1% SDS, 1 mM NaV, 1 mM NaF, en 1 EDTA-vrije protease inhibitor tablet (Roche)/10 mL in H2O) voor sonicatie met een Bioruptor Pico (Diagenode) gedurende 5 min (30 s op, 30 s uit). Lysaat werd gecentrifugeerd bij maximale snelheid gedurende 5 min om puin te verwijderen. Negen volumes van IP Verdunning Buffer (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, en 1 EDTA-vrije protease inhibitor tablet (Roche)/10 mL in H2O) werden toegevoegd aan het supernatant. In totaal werd 50 μl Dynabeads Protein G toegevoegd en werd het monster bij 4 °C gedurende 30 min. al roterend geïncubeerd. Eén procent van het monster werd als input bewaard en het resterende monster werd verdeeld over 3 buisjes. In totaal werd 50 μl Dynabeads Protein G geconjugeerd met 15 μl van het juiste antilichaam toegevoegd aan elke buis vóór overnachting incubatie bij 4 ° C, roterend. De met korrels gebonden complexen werden telkens 5 min. gewassen in 1 ml zoutarme buffer, zoutrijke buffer, LiCl-buffer en tweemaal met TE.
Om de met korrels gebonden complexen te elueren, werden de korrels geresuspendeerd in 50 μL elutiebuffer (100 mM NaHCO3, 1% w/v SDS) en gedurende 15 min. bij 65 °C geïncubeerd onder schudden bij 1000 omwentelingen per minuut op een thermomixer (Thermo Scientific). Elutie werd herhaald een tweede keer, en vervolgens 100 pi RNase Buffer (12 pi van 5 M NaCl, 0,2 pi 30 mg / ml RNase, en 88 pi TE) werd toegevoegd aan elke ChIP en input monster. De monsters werden geïncubeerd bij 37 ° C gedurende 20 min, gevolgd door de toevoeging van 100 pi van proteinase K buffer (2,5 pi 20 mg / ml proteinase K, 5 pi 20% SDS, en 92,5 pi TE) overnacht bij 65 ° C. Een gelijk volume van fenol:chloroform oplossing werd toegevoegd en grondig gemengd. Het mengsel werd overgebracht in MaXtract High Density-buisjes (Qiagen) en gedurende 8 min. bij 15 000 rpm gecentrifugeerd. De bovenste fase werd overgebracht naar nieuwe buisjes en gemengd met 1,5 pi 20 mg / ml glycogeen, 30 pi 3M natriumacetaat, en 800 pi ethanol. De monsters werden geïncubeerd bij – 80 °C tot bevroren en vervolgens gecentrifugeerd bij 15.000 rpm gedurende 30 min bij 4 °C. Het supernatant werd verwijderd en pellets werden gewassen in 800 pi 70% ijskoude ethanol en 10 min. gesponnen bij 4 ° C bij 15.000 rpm. Na zorgvuldige verwijdering van ethanol, werden pellets aan de lucht gedroogd en geresuspendeerd in 30 pi 10 mM Tris bij pH 8.
IP en input-DNA werden vervolgens gekwantificeerd met behulp van een Qubit 3.0 fluorometer. Bibliotheken werden bereid met behulp van de KAPA HyperPrep Kit (KK8505) en gesequenced met een Illumina NextSeq 500 tot een gemiddelde diepte van 28 miljoen reads per sample.
RNA-seq profilering
RNA werd geïsoleerd uit 3 miljoen cellen per monster met behulp van de Bio-Rad Aurum ™ Total RNA Mini Kit en gekwantificeerd met de Agilent RNA 6000 Nano Kit met de Agilent Bioanalyzer. Bibliotheken werden bereid door rRNA depletie met behulp van de Illumina TruSeq ® Stranded mRNA Library Prep Kit voor een lage concentratie van het uitgangsmonster en gesequeneerd door single-end sequencing op een Illumina NextSeq 500 tot een gemiddelde diepte van 18 miljoen leest per sample.
DNA methylering profilering
Genomisch DNA werd geïsoleerd met behulp van de AllPrep DNA/RNA Micro Kit (Qiagen). Om genoom-brede DNA-methylering status te beoordelen, voerden we mRRBS . Na fluorometrische kwantificering met behulp van een Qubit 3.0-instrument, we gedigesteerd genomisch DNA met het restrictie-enzym MspI (New England Biolabs) en grootte geselecteerd voor fragmenten ongeveer 100-250 basenparen in lengte met behulp van vaste fase omkeerbare immobilisatie (SPRI) kralen (MagBio Genomics). Het resulterende DNA onderging bisulfietconversie met behulp van de EZ DNA Methylation-Lightning Kit (Zymo Research). We creëerden bibliotheken van bisulfiet omgezet enkelstrengs DNA met behulp van de Pico Methyl-Seq Library Prep Kit (Zymo Research), die vervolgens werden gepoold voor sequencing op een Illumina NextSeq 500-instrument met behulp van de NextSeq 500/550 V2 High Output reagens kit (1 × 75 cycli) tot een minimale leesdiepte van 50 miljoen leest per monster.
Whole genome sequencing
Drie miljoen cellen van cellijnen of patiëntenmonsters werden gepelleteerd en geresuspendeerd in 1 ml Cell Lysis Solution (Qiagen) gemengd met 500 μg RNase A. De lysisreactie werd gedurende 15 min. bij 37 °C uitgevoerd. In totaal werd 333 pi eiwitprecipitatieoplossing (Qiagen) toegevoegd aan elk monster, dat vervolgens werd gedraaid en vervolgens gecentrifugeerd bij 2000 × g gedurende 10 minuten. Het supernatant werd gemengd met 1 mL isopropanol tot DNA-strengen neergeslagen uit de oplossing. Na het weggooien van het supernatant werd de DNA-pellet gewassen met 1 mL 70% ethanol en gecentrifugeerd bij 2000×g gedurende 1 min. De ethanol werd vervolgens uitgegoten en de pellet werd aan de lucht gedroogd gedurende 15 min voor resuspensie in 50 tot 100 ul van DNA-hydratatie-oplossing (Qiagen). DNA werd gesequeneerd met gepaarde-end Illumina sequencing bij 30× dekking.
Immunoprecipitatie
Een totaal van 100 miljoen cellen voor elke immunoprecipitatiereactie werden gepelleteerd en geïncubeerd in Buffer A (10 mM HEPES pH 8.0, 1,5 mM MgCl2, 10 mM KCl, 0,5 mM DTT) gedurende 10 min op ijs. Cellen werden vervolgens gelyseerd na 12 slagen met een 7-mL losse stamper weefselmolen (Wheaton, 357542) en gecentrifugeerd bij 2000 rpm gedurende 7 min. Nucleaire pellets werden geresuspendeerd in 5 volumes TENT buffer (50 mM Tris pH 7,5, 5 mM EDTA, 150 mM NaCl, 1% Triton X-100, 5 mM MgCl2) en behandeld met benzonase gedurende 30 min voor 5 passages door een 25 g × 5/8 in. spuit. De onoplosbare fractie werd verwijderd na centrifugatie bij 2000 rpm gedurende 7 min en ’s nachts geïncubeerd met Dynabeads Protein G gehybridiseerd met antilichaam. In totaal werden 2 miljoen cellen verwijderd voor invoer. De korrels en kernlysaten werden 6 keer gewassen met TENT-buffer en vervolgens geëlueerd in 0,1 M glycine pH 2,5 met 100 mM Tris pH 8,0 vooraf. NuPAGE LDS monsterbuffer werd toegevoegd aan eluaten en ingangen, die vervolgens werden geïncubeerd bij 70 ° C gedurende 15 min vóór analyse door western blot.
Public data collection
Public CTCF ChIP-seq gegevens werden verzameld uit Cistrome Data Browser (voor piek bestanden) en NCBI GEO (voor fastq bestanden, Additional file 2: tabel S1). Histon-modificatie ChIP-seq gegevens werden verzameld van NCBI GEO en ENCODE (voor bam-bestanden). Openbare RNA-seq gegevens in meerdere celtypes werden verzameld van ENCODE (voor fastq bestanden). DNA-methylering profilering gegevens werden verzameld van ENCODE (voor bedMethyl bestanden) en NCBI GEO. Hi-C-gegevens werden verzameld bij NCBI GEO en ENCODE (voor fastq-bestanden). ATAC-seq gegevens werden verzameld van NCBI GEO (voor fastq bestanden). Whole genome sequencing gegevens voor BRCA, COAD, LUAD, en PRAD monsters werden verzameld van International Cancer Genome Consortium (ICGC) Data Portal . Gedetailleerde informatie, waaronder de toetreding ID’s van alle openbare datasets verzameld in dit werk is te vinden in Additional file 6: Tabel S5.
Data processing
ChIP-seq data-analyse
Sequentie-uitlijning voor ChIP-seq gegevens in fastq bestanden werd uitgevoerd met behulp van dezelfde standaard analyse pijplijn zoals gebruikt in Cistrome DB , voor consistentie en reproduceerbaarheid. Alle sequentie data genomische uitlijning werden uitgevoerd met behulp van de Chilin pijplijn met standaard parameters ($ chilin eenvoudig -p smal -s hg38 –threads 8 -t IN.fq -i PRENAME -o OUTDIR). Kort samengevat werden de sequentieresets uitgelijnd aan het humane referentiegenoom (GRCH38/hg38) met behulp van BWA ($ 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). Sam-bestanden werden vervolgens omgezet in bam-bestanden met behulp van samtools ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). Voor CTCF ChIP-seq datasets werd MACS2 gebruikt om pieken op te roepen onder de FDR-drempel van 0,01 ($ macs2 callpeak –SPMR -B -q 0,01 –keep-dup 1 -g hs -t PRENAME.bam -n PRENAME –outidr OUTDIR). Pieken met een vouw verrijking van ten minste 4 werden behouden. Bigwiggle-bestanden werden gegenereerd met behulp van BEDTools en UCSC tools ($ 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). Ten slotte werden alleen de CTCF ChIP-seq monsters die ten minste 2000 pieken werden opgenomen in de downstream-integratieve analyse.
ATAC-seq data-analyse
Trim Galore werd gebruikt om de ruwe sequencing leest trimmen ($ trim_galore –nextera –phred33 –fastqc –paired R1.fq R2.fq -o OUTDIR). Leest werden uitgelijnd op het menselijk referentiegenoom (GRCH38/hg38) met behulp van Bowtie2 ($ bowtie2 -p 10 -X 2000 -x INDEX -1 R1.fq -2 R2.fq -S PRENAME.sam). Sam-bestanden werden vervolgens omgezet in bam-bestanden met behulp van samtools ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). Bedtools werd gebruikt om bam bestanden om te zetten in bed formaat ($ bamToBed -i PRENAME.bam -bedpe > PRENAME_PE.bed). Leest toegewezen aan mitochondria DNA werden verwijderd uit downstream-analyse.
RNA-seq data-analyse
RNA-seq datasets werden verwerkt met behulp van Salmon ($ salmon quant –gcBias -i INDEX -l A -p 8 {-1 R1.fq -2 R2.fq| -r IN.fq} -o OUTDIR). Transcriptoom index werd gebouwd op het menselijk referentiegenoom (GRCH38/hg38). Transcript-niveau abundantie schattingen werden samengevat tot het gen niveau met behulp van de “tximport” pakket voor differentiële expressie analyse. DESeq2 werd gebruikt om differentieel tot expressie komende genen te identificeren, en verschillende drempels gebruikt in verschillende analyses werden dienovereenkomstig vermeld in het manuscript.
Hi-C data-analyse
Hi-C gegevens werden verwerkt met HiC-Pro ($ HiC-Pro -i INDIR -o OUTDIR -c CONFIG -p). Contact kaarten werden gegenereerd met een resolutie van 5 kb. Ruwe matrixgegevens werden genormaliseerd met behulp van de aanpak beschreven in Normalization of Chromatin Interactions.
DNA methylatie data-analyse
DNA methylatie gegevens (voor T-ALL cellijnen en T-ALL patiënten) werden gedemultiplexed met bcl2fastq gevolgd door trimmen van 10 basisparen van de 5′ einde naar primer en adaptor sequenties met TrimGalore te verwijderen. Sequence alignment aan de GRCh38 / Hg38 referentie-genoom en methylatie oproepen werden uitgevoerd met Bismark ($ bismark –multicore 8 –bowtie2 -q -N 1 INDEX INFILE.fq). Dekking (tellingen) bestanden voor cytosines in CpG context werden gegenereerd met Bismark ($ bismark_methylation_extractor –multicore 8 –comprehensive –bedGraph INFILE_bismark_bt2.bam).
Whole genome sequencing data-analyse
Mutaties werden geïdentificeerd voor twee T-ALL cellijnen (Jurkat en CUTLL1) en twee T-ALL patiëntmonsters van de whole genome sequencing gegevens. We hebben de Illumina short-read sequenties uitgelijnd met het menselijke referentie genoom (GRCH38/hg38) met behulp van BWA mem. We gebruikten SAMBlaster om de discordante paren te identificeren, splitsen leest, en markeren de vermeende PCR duplicaten. We gebruikten SAMBAMBA om de SAM aligned om te zetten in het BAM formaat, en samtools werd gebruikt om de aligned te sorteren om een BAM bestand te creëren dat overeenkomt met elk monster.
We gebruikten VarDict om de varianten te identificeren die de unie CTCF bindingsplaatsen overlapten. We gebruikten alle standaard parameters, behalve “-f 0,1”, die werd gebruikt om varianten die werden ondersteund door meer dan 10% van de leest op die locatie te identificeren. We annoteerden de varianten met behulp van Variant Effect Predictor (VEP) en gebruikten aangepaste scripts om de varianten die TF binding beïnvloeden identificeren.
We gebruikten opnieuw VarDict om de varianten in de CTCF en NOTCH1 genen voor de vier monsters te identificeren. We gebruikten alle standaard parameters behalve “-f 0.1”, die werd gebruikt om varianten te identificeren die werden ondersteund door meer dan 10% van de gelezen op die locatie. We annoteerden de varianten met behulp van Variant Effect Predictor (VEP) , en vervolgens gefilterd om de mutaties die ofwel (a) niet gezien in meer dan 1% van een normale menselijke bevolking, of (b) had een CADD score van deleteriousness > 20, of (c) aanwezig was in de COSMIC database te identificeren.
Integratieve modellering en statistische analyse
Identificatie van CTCF-bindend repertoire in het menselijk genoom
Voor CTCF ChIP-seq verzamelden we in totaal 793 datasets, waaronder 787 openbare datasets en 6 datasets die we genereerden (Additional file 2: Tabel S1). In totaal werden 771 CTCF ChIP-seq datasets met pieken meer dan 2000 gebruikt in deze studie. Elke dataset kan MACS2-geïdentificeerde CTCF pieken opleveren in het bereik tussen 2050 en 198.021, met een mediaan van 46.451 en een totaal van 36.873.077 pieken (Additional file 1: Fig. S1a). De verdeling van de intervallengtes tussen aangrenzende CTCF piek toppen van alle 36.873.077 pieken uit de 771 datasets heeft een buigpunt bij ~ 150 bp (Additional file 1: Fig. S1c) die de grens aangeeft tussen dezelfde bindingsplaats en verschillende bindingsplaatsen . Daarom gebruikten we 150 bps als de cutoff om CTCF pieken samen te voegen. In de praktijk breidden we ± 75 bps uit van elke piektop om een 150 bp regio te genereren gecentreerd op de top om elke piek te vertegenwoordigen en voegden alle overlappende piekregio’s samen om een unie set van CTCF bindende sites te genereren, die 688.429 niet-overlappende sites bevat. Aan elke bindingsplaats werd een CTCF bezettingsscore toegekend, gedefinieerd als het totaal van ChIP-seq datasets die een piek vertonen binnen de plaats. Dienovereenkomstig definieerden we de bezettingsfrequentie als de verhouding van de bezettingsscore over het totale aantal CTCF ChIP-seq datasets. Om de robuustheid van de geïdentificeerde CTCF bindingsplaatsen verder te verzekeren, selecteerden we 285.467 hoog-vertrouwde plaatsen met bezettingsscore ≥ 3 voor downstream analyses. CTCF motieven binnen de unie bindingsplaatsen werden doorzocht door FIMO met Jaspar matrix (ID: MA0139.1), met een p-waardedrempel van 1e-4. Eén motief met de kleinste p-waarde werd behouden voor elke CTCF-bindingsplaats.
Identificatie van constitutieve CTCF-bindingsplaatsen
De verdeling van de bezettingsscores van alle 285.467 CTCF-bindingsplaatsen (Additional file 1: Fig. S1d, blauwe curve) laat zien dat de meerderheid van de CTCF-bindingsplaatsen in slechts enkele datasets voorkomt, en dat het aantal bindingsplaatsen afneemt met toenemende bezettingsscore wanneer de bezettingsscore klein is. Er zijn echter CTCF-bindingsplaatsen die zeer geconserveerd zijn in bijna alle datasets (b.v. bindingsplaatsen met een bezettingsscore van meer dan 600). We gebruiken een machtswet functie om de verdelingscurve (blauw) te passen die in Additional file 1: Fig. S1d om de cutoff voor constitutieve CTCF sites te bepalen. We duiden Oi als het aantal waargenomen CTCF bindende sites met bezettingsscore gelijk aan i, en Ei als het aantal verwachte CTCF sites met bezettingsscore gelijk aan i. De power law aanpassing aan gegevens Oi kan worden beschreven als (Additional file 1: Fig. S1d, groen):
We definiëren de cutoff A voor constitutieve CTCF-bindingsplaatsen als:
Met andere woorden, het totaal aan waargenomen CTCF-locaties met een bezettingsscore groter dan A zou 6 keer meer moeten zijn dan verwacht. We bepaalden vervolgens A = 615, en gebruikten een bezettingsfrequentie cutoff van 80% om 22.097 constitutieve CTCF-bindingsplaatsen te definiëren, wat overeenkomt met de bezettingsscore ≥ 616 in alle 771 CTCF ChIP-seq datasets.
Identificatie van kankerspecifieke gewonnen/verloren CTCF-bindingsplaatsen
We gebruikten de volgende 2 criteria om kankerspecifieke verloren CTCF-bindingsplaatsen te identificeren: (1) De CTCF-bindingsplaats moet een lagere bezettingsfrequentie hebben voor datasets van dat kankertype in vergelijking met de bezettingsfrequentie voor alle datasets en (2) het CTCF-bindingsniveau (gekwantificeerd als genormaliseerde ChIP-seq leestellingen) op de plaats is lager in kankerdatasets dan in andere datasets. Voor gewonnen CTCF sites gebruikten we de omgekeerde set criteria. Kort gezegd, voor elke CTCF bindingsplaats in elk type kanker, werd de bezettingsscore in de kanker datasets berekend samen met zijn bezettingsscore in alle 771 datasets. CTCF bindingsniveaus werden verkregen uit een genormaliseerde read count matrix waarin de ChIP-seq read counts (RPKM) eerst werden berekend voor unie CTCF bindingsplaatsen in alle datasets en daarna door kwantiel normalisatie. We gebruikten een ongepaarde tweestaart Student’s t test om het verschil in bindingsniveaus tussen verschillende groepen datasets te kwantificeren, en de p-waarde werd vervolgens aangepast met de Benjamini-Hochberg procedure. Bovendien werden de bindingsbezettingsscores en bindingsniveaus vergeleken tussen kankerdatasets en datasets van de gematchte normale weefsels of celtypes, om rekening te houden met de potentiële verstorende factor van weefselspecificiteit in plaats van kankerspecificiteit. Gedetailleerde criteria voor het identificeren van kankerspecifieke CTCF-bindingsplaatsen worden hieronder beschreven:
-
Kankerspecifieke verloren CTCF-bindingsplaatsen: (1) bezettingsfrequentie ≤ 0,2 in kankerdatasets; (2) bezettingsfrequentie ≥ 0,7 in 771 datasets; (3) bezettingsfrequentie ≥ 0.5 (met bezettingsscore ≥ 2) in gematchte normaalweefsel-datasets; (4) CTCF-niveaus zijn lager in kanker vergeleken met alle andere datasets (statistische score < 0), (5) CTCF-niveaus zijn lager in kanker vergeleken met gematchte normaalweefsel-datasets (statistische score < 0), (6) gemiddelde CTCF-bindingssignalen (RPKM) < 5 in kanker-datasets.
-
Kankerspecifieke gewonnen CTCF-bindingsplaatsen: (1) bezettingsfrequentie ≥ 0,5 (met bezettingsscore ≥ 2) in kankerdatasets, (2) bezettingsfrequentie ≤ 0,2 in 771 datasets, (3) bezettingsscore = 0 in gematchte normaalweefseldatasets, (4) CTCF-niveaus zijn significant hoger in kanker vergeleken met alle andere datasets (FDR ≤ 0.01), (5) CTCF-bindingsniveaus zijn significant hoger bij kanker in vergelijking met gematchte normaalweefsel-datasets (FDR ≤ 0,01), (6) gemiddelde CTCF-bindsignalen (RPKM) > 2 in kanker-datasets.
De specifieke gewonnen en verloren CTCF-bindingsplaatsen voor elk type kanker worden weergegeven in Extra bestand 4: Tabel S3.
Kwantificering van differentiële chromatine toegankelijkheid
We gebruikten de verwerkte gegevens van Ref. die een matrix van genormaliseerde ATAC-seq insertie telt binnen de TCGA pan-kanker piek set om de differentiële chromatine toegankelijkheid rond CTCF-bindingsplaatsen te beoordelen. Voor elk type kanker onder BRCA, CRC, LUAD, en PRAD, de pan-kanker ATAC-seq pieken die overlappen met geïdentificeerde kanker-type-specifieke verloren of gewonnen CTCF bindingsplaatsen werden gebruikt voor downstream analyses. De ATAC-seq differentiële score voor elke piek werd gekwantificeerd als de vouwverandering van het gemiddelde van de genormaliseerde ATAC-seq insertie tellingen van patiëntmonsters in het corresponderende kankertype versus van patiënten in andere kankertypes, en de ATAC-seq differentiële score werd vervolgens toegewezen aan de piek overlapte CTCF-bindingsplaats.
Voor consistentie, pasten we dezelfde aanpak toe die werd gebruikt voor TCGA ATAC-seq gegevens om de verzamelde ATAC-seq gegevens van T-ALL cellijn Jurkat en normale CD4 + T-cellen te analyseren. Een data matrix werd gegenereerd met behulp van ATAC-seq ruwe leestellingen op unie CTCF bindingsplaatsen voor alle Jurkat en T cel datasets. Kwantiel normalisatie werd toegepast op de log2 geschaald matrix (pseudo telling = 5). De ATAC-seq differentiële score werd gemeten als de vouw verandering van de gemiddelde genormaliseerde ATAC-seq telt tussen datasets van Jurkat versus CD4 + T-cel op elke CTCF-bindingsplaats.
Normalisatie van chromatine interacties
Gegeven een Hi-C contact kaart A = {aij}, de score aij weerspiegelt in kaart gebrachte leest tussen twee genomische regio’s i en j. Stel dat de bin grootte is 5 kb, regio’s i en j zal een genomische afstand van ∣i – j ∣ × 5kb hebben. Aangezien het contact waarschijnlijkheid tussen twee bakken afneemt met toenemende genomische afstand , hebben we genormaliseerd het contact kaart als volgt: voor een gegeven genomische afstand dk = k × 5kb, kwantificeren we een normalisatiefactor \( {overline{S}}_{d_k} \) als het gemiddelde interacties tussen alle bin paren met dezelfde genomische afstand dk in een hetzelfde chromosoom, bijv, \( {overline{S}_{d_k}=({\sum}_{j-i=k}{a}_{ij}right)/n \), waarbij n het totale aantal bin-paren met afstand dk is. De interactiescore aij tussen twee bins met afstand dk werd dan genormaliseerd door \( {{overline{S}}_{d_k} \) als \( {a}_{ij}^{\prime }={a}_{ij}/{{overline{S}}_{d_k} \). Met behulp van deze aanpak hebben we de matrix A genormaliseerd in A^{\prime }={a}_{ij}^{prime}right} \) binnen elk chromosoom.
Opsporing van differentiële chromatine interacties
We noemden de genormaliseerde Hi-C contact kaarten in de kanker dataset en de normale dataset als C = {cij} en N = {nij}, respectievelijk. Voor een bepaalde CTCF-bindingsplaats x (met coördinaat xc) en een vooraf gedefinieerde genomische afstand L, de chromatine interacties tussen x en de nabijgelegen niet-overlappende 5-kb bins met genomische afstand tot L worden verzameld uit C en N respectievelijk. Specifiek, interactie scores tussen x en de nabijgelegen 5-kb bakken in C worden verzameld als IC = {cij} terwijl ofwel i of j gelijk is aan ⌊xc/5kb⌋, en 0 < (j – i) × 5kb ≤ L. Evenzo werden de interactiescores tussen x en zijn nabije 5-kb bins in N verzameld als IN = {nij}. Een gepaarde two-tailed Student’s t-test werd vervolgens toegepast op IC en IN om de differentiële interactie tussen kanker en normale cellen rond CTCF-bindingsplaats x kwantificeren.
Binding van CTCF met genexpressie
In totaal werden 54 celtypen waarvoor zowel CTCF ChIP-seq gegevens en RNA-seq gegevens publiekelijk beschikbaar zijn geselecteerd (Additional file 6: Tabel S5) voor het onderzoeken van de associatie tussen CTCF binding en genexpressie voor elke CTCF-gen paar in hetzelfde chromosoom. Om het CTCF-bindende niveau te verkrijgen, werd een read count matrix gegenereerd met behulp van reads per kilobase per miljoen (RPKM) op unie CTCF bindende sites van ChIP-seq gegevens. De read count matrix werd geschaald met de vierkantswortel van RPKM gevolgd door kwantiel normalisatie. Genexpressie niveau werd gemeten voor elk gen met behulp van de vierkantswortel van transcripten per miljoen (TPM) van RNA-seq gegevens. Voor elk CTCF-gen-paar kwantificeerden we de associatie tussen de CTCF-locatie en het gen over alle 54 celtypes aan de hand van de correlatiecoëfficiënt R tussen het genexpressie-niveau en het genbindingsniveau (Fig. 3a). CTCF-gen paren werden geacht “sterk gecorreleerd” te zijn met een R2 groter dan 0,25, b.v, correlatiecoëfficiënt groter dan 0,5 of minder dan – 0,5, en de sterk gecorreleerde CTCF-genparen maken 1,3% uit van alle CTCF-genparen (Additional file 1: Fig. S8a).
Identificatie van constitutieve CTCF-gebonden chromatine domeinen
Voor elke CTCF-bindingsplaats, definieerden we de bijbehorende chromatine domein als de genomische regio die (1) omvat deze specifieke CTCF-bindingsplaats, (2) wordt begrensd door een paar constitutieve CTCF-bindingsplaatsen met motieven van tegengestelde oriëntaties, en (3) beslaat een minimum van 100 kb en een maximum van 1 MB regio aan elke kant van de CTCF-bindingsplaats. Figuur 3b bevat schematische weergave van hoe constitutieve CTCF-gebonden chromatine domeinen werden gedefinieerd.
Opsporing van DNA-methyleringsveranderingen rond CTCF-bindingsplaatsen
DNA-methyleringsveranderingen werden gedetecteerd binnen een 300-bp regio gecentreerd op elke CTCF-bindingsplaats. Regio’s met ten minste 3 CpGs gedekt door ten minste 5 gelezen (≥ 5×) in zowel kanker cellijnen en overeenkomstige normale weefsels werden behouden. Een 300-bp regio werd gedetecteerd als differentieel gemethyleerd als het gemiddelde differentiële methyleringsniveau van alle CpG’s (≥ 5×) binnen deze regio groter was dan 20%.
Opsporing van mutatiesnelheid en differentiële motief score
Voor elke CTCF-bindingsplaats werd het ruwe mutatiegetal berekend als het optreden van mutaties in alle monsters/patiënten bij elk afzonderlijk basenpaar binnen een 400-bp regio gecentreerd rond de CTCF-bindingsplaats. Het mutatiecijfer voor een groep CTCF-bindingsplaatsen werd berekend als het gemiddelde mutatiecijfer over het aantal CTCF-bindingsplaatsen voor elk basenpaar binnen de regio van 400 bp.
Motiefscore werd gemeten door de CTCF-positiegewichtsmatrix (Jaspar , Matrix ID: MA0139.1) aan een 19 bp DNA-sequentie gecentreerd rond het CTCF-motief of de CTCF-bindingsplaats te scoren met gebruikmaking van log-likelihoodratio’s (met achtergrondnucleotidenfrequentie als voor A,C,G,T). De differentiële motiefscore werd berekend door vergelijking van de motiefscores voor de referentie- en de gemuteerde sequenties.
DNA-sequentie motiefanalyse
DNA-sequentie motiefverrijkingsanalyse werd uitgevoerd met behulp van MDSeqPos (versie 1.0.0) op Cistrome met standaardparameters (-cisrome -Homo Sapien of Mus musculus). De novo motief analyses werden uitgevoerd met behulp van HOMER (versie 4.10) met findmotifs.pl module en MEME (versie 5.1.1) met de volgende parameters: meme -dna -mod zoops -maxw 20 -evt -0.01.
Identificatie van CTCF intra-domein differentieel geïnterageerde regio’s
Voor een bepaalde set van CTCF-bindingsplaatsen, werden de chromatine interactie veranderingen tussen een CTCF site en elk van zijn intra-domein niet-overlappende bins, gemeten uit genormaliseerde Hi-C contact kaarten in kankercellen over gematchte normale cellen, verzameld voor elk van de CTCF-bindingsplaatsen (Additional file 1: Fig. S14b). Regio’s met verminderde interacties (log2 FC < -1, gemiddelde log2 interactie > 0) met kanker-specifieke verloren CTCF-bindingsplaatsen, en regio’s met verhoogde interacties (log2 FC > 1, gemiddelde log2 interactie > 0) met kanker-specifieke gewonnen CTCF-bindingsplaatsen werden gebruikt voor stroomafwaartse transcriptiefactor (TF) verrijkingsanalyse.
Transcriptiefactor verrijkingsanalyse
Een herziene versie van het BART-algoritme werd gebruikt voor TF verrijkingsanalyse. Kort gezegd, een verzameling van unie DNase I overgevoelige plaatsen (UDHS) werd eerder gecureerd als een repertoire van alle kandidaat cis-regulerende elementen in het menselijk genoom, en 7032 ChIP-seq datasets werden verzameld voor 883 TFs , met elke TF met een of meer ChIP-seq datasets van meerdere celtypen of omstandigheden. Voor elke TF op UDHS werd een binair profiel gegenereerd dat aangeeft of de TF ten minste één piek van een van zijn ChIP-seq datasets heeft die zich binnen elk van de UDHS bevindt. Bindingsverrijkingsanalyse werd voor elke TF toegepast door de TF-binding op een subset van UDHS die de geselecteerde genomische regio’s overlappen, te vergelijken met de TF-binding op UDHS. p-waarde werd verkregen met behulp van Fisher’s exacte test met twee steekproeven.