Forsøgsprocedure
Patient xenografting og cellekultur
De humane T-ALL-cellelinjer omfatter CUTLL1 (gave fra Adolfo Ferrando, Columbia University) og JURKAT (American Type Culture Collection (ATCC), Manassas, VA, #CCL-119) . Cellerne blev dyrket i RPMI1640-medium med l-glutamin og 25 mM HEPES (Corning) suppleret med 10% varmeinaktiveret føtalt kvægserum (Sigma-Aldrich), 10 U/mL penicillin-streptomycin (Gibco) og 1× glutaMAX (Gibco) i en befugtet inkubator ved 37 °C og 5% CO2. Cellerne testes regelmæssigt for tilstedeværelse af mycoplasmaer ved hjælp af Lonza Walkersville MycoAlert Mycoplasma Detection Kit (sidste test i januar 2020). Cellelinjerne holdes i kultur i højst 20 passager og autentificeres ved hjælp af kort-tandem gentagelsesprofilering (JURKAT) eller ved hjælp af PCR til påvisning af TCRb-NOTCH1-translokationen (TCRBJ2S4CUTLL1F:5′-GGACCCGGCTCTCTCAGAGTGCT-3′, NOTCH1CUTTL1R:5′-TCCCGCCCTCCACAAAATAAGG-3′). Den sidste celleautentificering blev foretaget i februar 2020. Humane CD4+ T-celler blev købt hos AllCells. Primære humane prøver blev indsamlet med informeret samtykke og analyseret under tilsyn af det institutionelle bedømmelsesudvalg ved universitetet i Padova, Associazone Italiana di Ematologia e Oncologia Pediatrica og de kliniske pædiatriske ALL 2000/2006-prøver i Berlin-Frankfurt-Münster (AIEOP-BFM). Der blev indhentet informeret samtykke til at bruge overskydende materiale til forskningsformål fra alle patienterne ved forsøgets start i overensstemmelse med Helsinki-erklæringen.
Antistoffer og reagenser
Vestern blots blev udført ved hjælp af følgende antistoffer: Der blev anvendt følgende antistofferne: Aktin og CTCF fra Millipore Sigma (klon C4; 07-729) og kløvet NOTCH1 (Val1744) fra Cell Signaling Technology (4147). ChIP-seq blev udført ved hjælp af følgende antistoffer: CTCF fra Millipore Sigma (07-729), H3K27Ac (8173S) og H3K27me3 (9733S) fra Cell Signaling Technology og H3K4me1 (07-473) fra Millipore.
In situ Hi-C
In situ Hi-C blev udført på CD4+ T-celler, Jurkat, CUTLL1 og patient xenografts som tidligere beskrevet . Kort fortalt blev cellerne krydsbundet med 1% formaldehyd i 10 min ved stuetemperatur. Pr. Hi-C-reaktion blev 5 millioner celler lyseret, og kerner blev permeabiliseret. DNA blev fordøjet med MboI fra New England Biolabs (R0147M). De fordøjede fragmenter blev mærket med biotinyleret d-ATP fra Jena Bioscience (NU-835-BIO14-S) og ligeret. Efter RNase-behandling og Proteinase K-behandling for at ophæve tværbindinger blev kerner soniceret ved hjælp af en Covaris E220 for at opnå en gennemsnitlig fragmentlængde på 400 bp. Streptavidin perler fra Thermo Fisher Scientific (65001) blev anvendt til at trække biotinmærkede fragmenter ned. Efter rensning og isolering af DNA blev de endelige biblioteker fremstillet ved hjælp af NEBNext® Ultra™ II DNA Library Prep Kit for Illumina® og sekventeret ved parvis sekventering med en læselængde på 150 bp på en Illumina HiSeq 2500 for at producere i gennemsnit 400 millioner læs pr. prøve.
ChIP-seq-profilering
CD4+ T-celler, Jurkat, CUTLL1 og patient xenografts blev tværbundet med 1 % formaldehyd og 1 % føtal bovin serum i PBS i 10 min. ved stuetemperatur. Reaktionen blev slukket med 0,2 M glycin ved stuetemperatur i 5 min. Cellerne blev derefter vasket med PBS og pelleteret.
For CTCF ChIPs blev immunoprecipitation udført baseret på en protokol beskrevet tidligere . En pellet indeholdende 50 millioner celler blev lyseret med 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) i 10 min ved 4 °C. Kerner blev pelleteret ved 1350×g i 7 min og resuspenderet i 10 mM Tris pH 8, 1 mM EDTA og 0,1 % SDS. Kromatin blev skåret med et Covaris E220-system til en gennemsnitlig fragmentlængde på 400 bp og drejet ved 15.000 rpm i 10 min. for at fjerne uopløseligt kromatin og affald. Supernatanten blev inkuberet med 20 μL Dynabeads Protein G i 30 min, inden perlerne blev kasseret. En procent af det samlede volumen blev gemt som input, og resten blev inkuberet med anti-CTCF-antistof natten over. I alt blev der tilsat 100 μL Dynabeads Protein G i 2 h. De bundne fragmenter blev vasket to gange med 1 mL lav saltbuffer (20 mM Tris-HCl pH 8,0, 150 mM NaCl, 2 mM EDTA, 1 % w/v Triton X-100 og 0.1% w/v SDS), én gang med høj saltbuffer (20 mM Tris-HCl pH 8,0, 500 mM NaCl, 2 mM EDTA, 1% w/v Triton X-100 og 0,1% w/v SDS), én gang med lithiumchloridbuffer (10 mM Tris-HCl pH 8.0, 250 mM LiCl, 1 mM EDTA, 1 % w/v NP-40 og 1 % w/v deoxycholsyre) og to gange med TE (10 mM Tris pH 8, 1 mM EDTA).
For histone ChIP’er blev cellerne lyseret i 375 μL kerneinkubationsbuffer (15 mM Tris pH 7.5, 60 mM KCl, 150 mM NaCl, 15 mM MgCl2, 1 mM CaCl2, 250 mM saccharose, 0,3% NP-40, 1 mM NaV, 1 mM NaF og 1 EDTA-fri proteaseinhibitortablet (Roche)/10 mL i H2O) i 10 min på is i 10 min. Kerner blev vasket én gang med digestbuffer (10 mM NaCl, 10 mM Tris pH 7,5, 3 mM MgCl2, 1 mM CaCl2, 1 mM NaV, 1 mM NaF og 1 EDTA-fri proteaseinhibitortablet (Roche)/10 mL i H2O) og resuspenderet i 57-μL Digestbuffer indeholdende 4,5 enheder MNase (USB) i 1 time ved 37 °C. MNaseaktiviteten blev slukket i 10 min på is ved tilsætning af EDTA til en slutkoncentration på 20 mM. Kerner blev pelleteret og resuspenderet i 300-μL kernelysisbuffer (50 mM Tris-HCl pH 8,0, 10 mM EDTA pH 8,0, 1 % SDS, 1 mM NaV, 1 mM NaF og 1 EDTA-fri proteaseinhibitortablet (Roche)/10 mL i H2O) før sonicering med en Bioruptor Pico (Diagenode) i 5 min (30 s tændt, 30 s slukket). Lysatet blev centrifugeret ved maksimal hastighed i 5 min. for at fjerne debris. Ni volumener IPfortyndingsbuffer (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 og 1 EDTA-fri proteaseinhibitortablet (Roche)/10 mL i H2O) blev tilsat til supernatanten. I alt blev der tilsat 50 μL Dynabeads Protein G, og prøven blev inkuberet ved 4 °C i 30 min. under rotation. En procent af prøven blev beholdt som input, og den resterende prøve blev fordelt på 3 rør. I alt blev der tilsat 50 μL Dynabeads Protein G konjugeret til 15 μL af det relevante antistof til hvert rør før inkubation natten over ved 4 °C under rotation. Perlebundne komplekser blev vasket i 5 min hver i 1 mL lav-saltbuffer, høj-saltbuffer, LiCl-buffer og to gange med TE.
For at eluere perlebundne komplekser blev perlerne resuspenderet i 50 μL elueringsbuffer (100 mM NaHCO3, 1 % w/v SDS) og inkuberet ved 65 °C i 15 min under omrystning ved 1000 RPM på en termomixer (Thermo Scientific). Elueringen blev gentaget en anden gang, og derefter blev 100 μL RNase Buffer (12 μL 5 M NaCl, 0,2 μL 30 mg/mL RNase og 88 μL TE) tilsat til hver ChIP- og inputprøve. Prøverne blev inkuberet ved 37 °C i 20 min, efterfulgt af tilsætning af 100 μL proteinase K-buffer (2,5 μL 20 mg/mL proteinase K, 5 μL 20 % SDS og 92,5 μL TE) natten over ved 65 °C. En lige stor mængde phenol:chloroform-opløsning blev tilsat og blandet grundigt. Blandingen blev overført til MaXtract High Density-rør (Qiagen) og centrifugeret i 8 minutter ved 15 000 omdrejninger pr. minut. Den øverste fase blev overført til nye rør og blandet med 1,5 μL 20 mg/mL glykogen, 30 μL 3M natriumacetat og 800 μL ethanol. Prøverne blev inkuberet ved – 80 °C, indtil de blev frosset, og derefter centrifugeret ved 15 000 rpm i 30 min. ved 4 °C. Supernatanten blev fjernet, og pellets blev vasket i 800 μL 70 % iskold ethanol og centrifugeret i 10 min. ved 4 °C ved 15 000 rpm. Efter omhyggelig fjernelse af ethanol blev pellets lufttørret og resuspenderet i 30 μL 10 mM Tris ved pH 8.
IP og input-DNA blev derefter kvantificeret ved hjælp af et Qubit 3.0-fluorometer. Biblioteker blev fremstillet ved hjælp af KAPA HyperPrep Kit (KK8505) og sekventeret med en Illumina NextSeq 500 til en gennemsnitlig dybde på 28 millioner læsninger pr. prøve.
RNA-seq-profilering
RNA blev isoleret fra 3 millioner celler pr. prøve ved hjælp af Bio-Rad Aurum™ Total RNA Mini Kit fra Bio-Rad Aurum™ og kvantificeret med Agilent RNA 6000 Nano Kit med Agilent Bioanalyzer. Biblioteker blev fremstillet ved rRNA-depletion ved hjælp af Illumina TruSeq® Stranded mRNA Library Prep Kit til en lav koncentration af startprøven og sekventeret ved single end sekventering på en Illumina NextSeq 500 til en gennemsnitlig dybde på 18 millioner reads pr. prøve.
DNA-methyleringsprofilering
Genomisk DNA blev isoleret ved hjælp af AllPrep DNA/RNA Micro Kit (Qiagen). For at vurdere DNA-methyleringsstatus i hele genomet udførte vi mRRBS . Efter fluorometrisk kvantificering ved hjælp af et Qubit 3.0-instrument fordøjede vi genomisk DNA med restriktionsenzymet MspI (New England Biolabs) og størrelsesselekterede for fragmenter med en længde på ca. 100-250 basepar ved hjælp af fastfase reversibel immobilisering (SPRI) perler (MagBio Genomics). Det resulterende DNA blev bisulfitkonverteret ved hjælp af EZ DNA Methylation-Lightning Kit (Zymo Research). Vi skabte biblioteker fra bisulfitkonverteret enkeltstrenget DNA ved hjælp af Pico Methyl-Seq Library Prep Kit (Zymo Research), som derefter blev samlet til sekventering på et Illumina NextSeq 500-instrument ved hjælp af NextSeq 500/550 V2 High Output-reagenskit (1 × 75 cyklusser) til en minimum læsedybde på 50 millioner læsninger pr. prøve.
Helgenomsekventering
Tre millioner celler fra cellelinjer eller patientprøver blev pelleteret og resuspenderet i 1 mL Cell Lysis Solution (Qiagen) blandet med 500 μg RNase A. Lysereaktionen blev udført ved 37 °C i 15 minutter. I alt blev der tilsat 333 μL proteinudfældningsopløsning (Qiagen) til hver prøve, som derefter blev vortexet og derefter centrifugeret ved 2000×g i 10 min. Supernatanten blev blandet med 1 mL isopropanol, indtil DNA-strenge udfældede fra opløsningen. Efter at have kasseret supernatanten blev DNA-pellet vasket med 1 mL 70 % ethanol og centrifugeret ved 2000×g i 1 min. Ethanolet blev derefter hældt ud, og pellet blev lufttørret i 15 min., inden det blev resuspenderet i 50-100 μL DNA-hydreringsopløsning (Qiagen). DNA blev sekventeret med parvis-end Illumina-sekventering ved 30× dækning.
Immunopræcipitation
I alt 100 millioner celler til hver immunopræcipitationsreaktion blev pelleteret og inkuberet i buffer A (10 mM HEPES pH 8,0, 1,5 mM MgCl2, 10 mM KCl, 0,5 mM DTT) i 10 min på is. Cellerne blev derefter lyseret ved 12 slag med en 7-mL vævsmølle med løs pestel (Wheaton, 357542) og centrifugeret ved 2000 rpm i 7 min. Kernepellets blev resuspenderet i 5 volumener TENT-buffer (50 mM Tris pH 7,5, 5 mM EDTA, 150 mM NaCl, 1% Triton X-100, 5 mM MgCl2) og behandlet med benzonase i 30 min. før 5 passager gennem en 25 g × 5/8 tommer sprøjte. Den uopløselige fraktion blev fjernet efter centrifugering ved 2000 rpm i 7 min og inkuberet natten over med Dynabeads Protein G hybridiseret med antistof. I alt 2 millioner celler blev fjernet til input. Beads og kerne-lysater blev vasket 6 gange med TENT-buffer og derefter elueret i 0,1 M glycin pH 2,5 med 100 mM Tris pH 8,0 forud. NuPAGE LDS-prøvebuffer blev tilsat eluater og input, som derefter blev inkuberet ved 70 °C i 15 min. før analyse ved western blot.
Offentlig dataindsamling
Offentlige CTCF ChIP-seq-data blev indsamlet fra Cistrome Data Browser (for peak-filer) og NCBI GEO (for fastq-filer, Additional file 2: Table S1). Histonmodifikations ChIP-seq-data blev indsamlet fra NCBI GEO og ENCODE (for bam-filer). Offentlige RNA-seq-data i flere celletyper blev indsamlet fra ENCODE (for fastq-filer). Data om DNA-methyleringsprofilering blev indsamlet fra ENCODE (for bedbedMethyl-filer) og NCBI GEO. Hi-C-data blev indsamlet fra NCBI GEO og ENCODE (for fastq-filer). ATAC-seq-data blev indsamlet fra NCBI GEO (for fastq-filer). Helgenom-sekventeringsdata for BRCA-, COAD-, LUAD- og PRAD-prøver blev indsamlet fra International Cancer Genome Consortium (ICGC) Data Portal . Detaljerede oplysninger, herunder accession ID’er for alle offentlige datasæt indsamlet i dette arbejde, kan findes i Additional file 6: Table S5.
Data processing
ChIP-seq dataanalyse
Sequence alignment for ChIP-seq data i fastq filer blev udført ved hjælp af den samme standard analyse pipeline som anvendt i Cistrome DB , for konsistens og reproducerbarhed. Alle sekvensdata genomisk alignment blev udført ved hjælp af Chilin-pipeline med standardparametre ($ chilin simple -p narrow -s hg38 –threads 8 -t IN.fq -i PRENAME -o OUTDIR). Kort fortalt blev sekvenslæsninger afstemt med det menneskelige referencegenom (GRCH38/hg38) ved hjælp af 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-filerne blev derefter konverteret til bam-filer ved hjælp af samtools ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). For CTCF ChIP-seq-datasæt blev MACS2 anvendt til at kalde toppe under FDR-tærsklen på 0,01 ($ macs2 callpeak –SPMR -B -q 0,01 –keep-dup 1 -g hs -t PRENAME.bam -n PRENAME –outidr OUTDIR). Peaks med en foldberigelse på mindst 4 blev bevaret. Bigwiggle-filer blev genereret ved hjælp af BEDTools og UCSC-værktøjer ($ 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). Endelig blev kun de CTCF ChIP-seq-prøver, der har mindst 2000 peaks, inkluderet i den efterfølgende integrative analyse.
ATAC-seq-dataanalyse
Trim Galore blev brugt til at trimme de rå sekventeringslæsninger ($ trim_galore –nextera –phred33 –fastqc –paired R1.fq R2.fq -o OUTDIR). Læserne blev afstemt med det menneskelige referencegenom (GRCH38/hg38) ved hjælp af Bowtie2 ($ bowtie2 -p 10 -X 2000 -x INDEX -1 R1.fq -2 R2.fq -S PRENAME.sam). Sam-filerne blev derefter konverteret til bam-filer ved hjælp af samtools ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). Bedtools blev brugt til at konvertere bam-filer til bed-format ($ bamToBed -i PRENAME.bam -bedpe > PRENAME_PE.bed). Læsninger, der var kortlagt til mitokondrie-DNA, blev udeladt fra downstream-analyse.
RNA-seq-dataanalyse
RNA-seq-datasæt blev behandlet ved hjælp af Salmon ($ salmon quant –gcBias -i INDEX -l A -p 8 {-1 R1.fq -2 R2.fq| -r IN.fq} -o OUTDIR). Transkriptomindekset blev bygget på det menneskelige referencegenom (GRCH38/hg38). Estimater af transkriptniveauets abundans blev sammenfattet til genniveauet ved hjælp af “tximport”-pakken til differentiel ekspressionsanalyse. DESeq2 blev anvendt til at identificere differentielt udtrykte gener, og forskellige tærskelværdier, der blev anvendt i forskellige analyser, blev opført tilsvarende i manuskriptet.
Hi-C-dataanalyse
Hi-C-data blev behandlet ved hjælp af HiC-Pro ($ HiC-Pro -i INDIR -o OUTDIR -c CONFIG -p). Der blev genereret kontaktkort med en opløsning på 5 kb. Rå matrixdata blev normaliseret ved hjælp af den fremgangsmåde, der er beskrevet i Normalisering af kromatininteraktioner.
DNA-methyleringsdataanalyse
DNA-methyleringsdata (for T-ALL cellelinjer og T-ALL patienter) blev demultiplexet med bcl2fastq efterfulgt af trimning af 10 basepar fra 5′-enden for at fjerne primer- og adaptor-sekvenser ved hjælp af TrimGalore . Sekvenstilpasning til GRCh38/hg38-referencegenomet og methyleringskald blev udført med Bismark ($ bismark –multicore 8 –bowtie2 -q -N 1 INDEX INFILE.fq). Dækningsfiler (counts) for cytosiner i CpG-kontekst blev genereret ved hjælp af Bismark ($ bismark_methylation_extractor –multicore 8 –comprehensive –bedGraph INFILE_bismark_bt2.bam).
Analyse af helgenom-sekventeringsdata
Mutationer blev identificeret for to T-ALL-cellelinjer (Jurkat og CUTLL1) og to T-ALL-patientprøver fra helgenom-sekventeringsdataene. Vi tilpassede Illumina short-read-sekvenserne til det menneskelige referencegenom (GRCH38/hg38) ved hjælp af BWA mem. Vi brugte SAMBlaster til at identificere de uoverensstemmende par, opdele læsninger og markere de formodede PCR-duplikater. Vi brugte SAMBAMBA til at konvertere SAM-justeret til BAM-formatet, og samtools blev brugt til at sortere de justerede for at skabe en BAM-fil svarende til hver prøve.
Vi brugte VarDict til at identificere de varianter, der overlappede union CTCF-bindingsstederne. Vi brugte alle standardparametre undtagen “-f 0,1”, som blev brugt til at identificere varianter, der blev understøttet af mere end 10 % af læsningerne på det pågældende sted. Vi annoterede varianterne ved hjælp af Variant Effect Predictor (VEP) og brugte brugerdefinerede scripts til at identificere de varianter, der påvirker TF-binding.
Vi brugte igen VarDict til at identificere varianterne i CTCF- og NOTCH1-generne for de fire prøver. Vi brugte alle standardparametre undtagen “-f 0,1”, som blev brugt til at identificere varianter, der blev understøttet af mere end 10 % af læsningerne på det pågældende sted. Vi annoterede varianterne ved hjælp af Variant Effect Predictor (VEP) , og filtrerede den derefter for at identificere de mutationer, der enten (a) ikke blev set i mere end 1 % af nogen normal menneskelig population, eller (b) havde en CADD-score for deleteriousness > 20, eller (c) var til stede i COSMIC-databasen.
Integrativ modellering og statistisk analyse
Identificering af CTCF-bindingsrepertoire i det menneskelige genom
For CTCF ChIP-seq indsamlede vi i alt 793 datasæt, herunder 787 offentlige datasæt og 6 datasæt, som vi genererede (Yderligere fil 2: Tabel S1). I alt blev 771 CTCF ChIP-seq-datasæt med toppe på mere end 2000 brugt i denne undersøgelse. Hvert datasæt kan give MACS2-identificerede CTCF-toppe i intervallet mellem 2050 og 198.021 med en median på 46.451 og i alt 36.873.077 toppe (Yderligere fil 1: Fig. S1a). Fordelingen af intervallængderne mellem tilstødende CTCF-toppe af alle 36.873.077 toppe fra de 771 datasæt har et bøjningspunkt ved ~ 150 bp (Yderligere fil 1: Fig. S1c), der angiver grænsen mellem det samme bindingssted og forskellige bindingssteder . Derfor brugte vi 150 bps som grænseværdi for sammenlægning af CTCF-toppe. I praksis udvidede vi ± 75 bps fra hver top for at generere en 150-bp-region centreret på toppen for at repræsentere hver top og fusionerede alle overlappende topregioner for at generere et unionssæt af CTCF-bindingssteder, som indeholder 688.429 ikke-overlappende steder. Hvert bindingssted blev tildelt en CTCF-belægningsscore, defineret som tallet af ChIP-seq-datasæt, der viser en top inden for stedet. I overensstemmelse hermed definerede vi belægningsfrekvensen som forholdet mellem belægningsscoren og det samlede antal CTCF ChIP-seq-datasæt. For yderligere at sikre robustheden af de identificerede CTCF-bindingssteder valgte vi 285.467 steder med høj tillid med belægningsscore ≥ 3 til downstream-analyser. CTCF-motiver inden for unionsbindingsstederne blev søgt ved hjælp af FIMO med Jaspar-matrix (ID: MA0139.1) med en p-værditærskel på 1e-4. Et motiv med den mindste p-værdi blev bevaret for hvert CTCF-bindingssted.
Identificering af konstitutive CTCF-bindingssteder
Distributionen af belægningsscore for alle 285.467 CTCF-bindingssteder (Yderligere fil 1: Fig. S1d, blå kurve) viser, at størstedelen af CTCF-bindingsstederne kun forekommer i få datasæt, og antallet af bindingssteder falder med stigende belægningsscore, når belægningsscoren er lille. Der er imidlertid CTCF-bindingssteder, som er meget bevaret på tværs af næsten alle datasæt (f.eks. bindingssteder med en belægningsscore på over 600). Vi bruger en power law-funktion til at tilpasse distributionskurven (blå), der er vist i Additional file 1: Fig. S1d for at bestemme cutoff for konstitutive CTCF-steder. Vi betegner Oi som antallet af observerede CTCF-bindingssteder med en belægningsscore lig med i, og Ei som antallet af forventede CTCF-steder med en belægningsscore lig med i. Effektlovstilpasningen til dataene Oi kan beskrives som (Additional file 1: Fig. S1d, grøn):
Vi definerer cutoff A for konstitutive CTCF-bindingssteder som:
Med andre ord bør det samlede antal observerede CTCF-steder med en belægningsscore større end A være 6 gange så stort som forventet. Vi bestemte derefter A = 615 og brugte en belægningsfrekvens cutoff på 80 % til at definere 22.097 konstitutive CTCF-bindingssteder, hvilket svarer til belægningsscore ≥ 616 i alle 771 CTCF ChIP-seq-datasæt.
Identificering af kræftspecifikke vundne/tabte CTCF-bindingssteder
Vi brugte følgende 2 kriterier til at identificere kræftspecifikke tabte CTCF-bindingssteder: (1) CTCF-bindingsstedet skal have en lavere belægningsfrekvens for datasæt af den pågældende kræfttype sammenlignet med belægningsfrekvensen for alle datasæt, og (2) CTCF-bindingsniveauet (kvantificeret som normaliserede ChIP-seq-læsningstællinger) på stedet er lavere i kræftdatasæt end i andre datasæt. For opnåede CTCF-steder anvendte vi det omvendte sæt af kriterier. Kort sagt blev for hvert CTCF-bindingssted i hver kræfttype belægningsscoren i kræftdatasættene beregnet sammen med dens belægningsscore i alle 771 datasæt. CTCF-bindingsniveauer blev opnået fra en normaliseret læsetællingsmatrix, hvor ChIP-seq-læsetællingerne (RPKM) først blev beregnet for fagforeningens CTCF-bindingssteder i alle datasæt og derefter efterfulgt af kvantilnormalisering. Vi brugte uparret tosidet Student’s t-test til at kvantificere forskellen i bindingsniveauer mellem forskellige grupper af datasæt, og p-værdien blev derefter justeret ved hjælp af Benjamini-Hochberg-proceduren . Desuden blev bindingsbelægningsscore og bindingsniveauer sammenlignet mellem kræftdatasæt og datasæt fra de matchede normale vævs- eller celletyper for at tage hensyn til den potentielle forvirrende faktor af vævsspecificitet snarere end kræftspecificitet. Detaljerede kriterier for identifikation af kræftspecifikke CTCF-bindingssteder er beskrevet nedenfor:
-
Kræftspecifikke tabte CTCF-bindingssteder: (1) belægningsfrekvens ≤ 0,2 i kræftdatasæt; (2) belægningsfrekvens ≥ 0,7 i 771-datasæt; (3) belægningsfrekvens ≥ 0.5 (med belægningsscore ≥ 2) i matchede datasæt af normalt væv; (4) CTCF-niveauer er lavere i kræft sammenlignet med alle andre datasæt (statistisk score < 0), (5) CTCF-niveauer er lavere i kræft sammenlignet med matchede datasæt af normalt væv (statistisk score < 0), (6) gennemsnitlige CTCF-bindingssignaler (RPKM) < 5 i kræftdatasæt.
-
Kræftspecifikke opnåede CTCF-bindingssteder: (1) belægningsfrekvens ≥ 0,5 (med belægningsscore ≥ 2) i kræftdatasæt, (2) belægningsfrekvens ≤ 0,2 i 771 datasæt, (3) belægningsscore = 0 i matchede datasæt af normalt væv, (4) CTCF-niveauer er signifikant højere i kræft sammenlignet med alle andre datasæt (FDR ≤ 0.01), (5) CTCF-bindingsniveauer er signifikant højere i kræft sammenlignet med matchede datasæt af normalt væv (FDR ≤ 0,01), (6) gennemsnitlige CTCF-bindingssignaler (RPKM) > 2 i kræftdatasæt.
De specifikke vundne og tabte CTCF-bindingssteder for hver kræfttype er vist i Additional file 4: Tabel S3.
Kvantificering af differentiel kromatintilgængelighed
Vi brugte de behandlede data fra Ref. der omfatter en matrix af normaliserede ATAC-seq-indføjningstællinger inden for TCGA pan-cancer peak-sættet til at vurdere den differentielle kromatintilgængelighed omkring CTCF-bindingssteder. For hver kræfttype blandt BRCA, CRC, LUAD og PRAD blev de pan-kræft ATAC-seq-toppe, der overlapper med identificerede kræfttypespecifikke tabte eller vundne CTCF-bindingssteder, anvendt til downstream-analyser. ATAC-seq differential score for hver top blev kvantificeret som foldændringen af gennemsnittet af de normaliserede ATAC-seq indsætningstællinger fra patientprøver i den tilsvarende kræfttype i forhold til patienter i andre kræfttyper, og ATAC-seq differential score blev derefter tildelt toppen overlappede CTCF-bindingssted.
For konsistensensens skyld anvendte vi den samme fremgangsmåde, der blev anvendt til TCGA ATAC-seq-data, til at analysere de indsamlede ATAC-seq-data fra T-ALL cellelinjen Jurkat og normale CD4+ T-celler. Der blev genereret en datamatrix ved hjælp af ATAC-seq rå læsningstællinger på union CTCF-bindingssteder for alle Jurkat- og T-celle-datasæt. Der blev anvendt kvantilenormalisering på den log2-skalerede matrix (pseudotælling = 5). ATAC-seq-differentialscore blev målt som foldændringen af de gennemsnitlige normaliserede ATAC-seq-tællinger mellem datasæt af Jurkat versus CD4+ T-celler på hvert CTCF-bindingssted.
Normalisering af kromatininteraktioner
Givet et Hi-C-kontaktkort A = {aij} afspejler scoren aij kortlagte læsninger mellem to genomiske regioner i og j. Antag, at bin-størrelsen er 5 kb, vil regionerne i og j have en genomisk afstand på ∣i – j ∣ × 5kb. Da kontaktsandsynligheden mellem to bin falder med stigende genomisk afstand , normaliserede vi kontaktkortet på følgende måde: For en given genomisk afstand dk = k × 5kb kvantificerer vi en normaliseringsfaktor \( { {\overline{S}}}_{d_k} \) som de gennemsnitlige interaktioner mellem alle bin-par med samme genomiske afstand dk i et og samme kromosom, f.eks, \( {\overline{S}}}_{d_k}=\left({\sum}_{j-i=k}{a}_{ij}\right)/n \), hvor n er det samlede antal bin-par med afstand dk. Interaktionsscoren aij mellem to bins med afstand dk blev derefter normaliseret med \( {\overline{S}}}_{d_k} \) som \( {a}_{ij}^{\prime }={a}_{ij}/{\overline{S}}_{d_k} \). Ved hjælp af denne fremgangsmåde normaliserede vi matrixen A til \( A^{\prime }=\left\{{a}_{ij}^{\prime}\right\} \) inden for hvert kromosom.
Detektion af differentielle kromatininteraktioner
Vi betegner de normaliserede Hi-C-kontaktkort i kræftdatasættet og det normale datasæt som henholdsvis C = {cij} og N = {nij}. For et givet CTCF-bindingssted x (med koordinat xc) og en foruddefineret genomisk afstand L indsamles kromatininteraktionerne mellem x og dets nærliggende ikke-overlappede 5-kb-bins med genomisk afstand op til L fra henholdsvis C og N. Interaktionsscorerne mellem x og de nærliggende 5-kb-bins i C indsamles specifikt som IC = {cij} , mens enten i eller j er lig med ⌊xc/5kb⌋, og 0 < (j – i) × 5kb ≤ L. På samme måde blev interaktionsscorerne mellem x og dets nærliggende 5-kb-bins i N indsamlet som IN = {nij}. Der blev derefter anvendt en parret tosidet Student’s t-test på IC og IN for at kvantificere den differentielle interaktion mellem kræftceller og normale celler omkring CTCF-bindingsstedet x.
Sammenhæng mellem CTCF-binding og genekspression
I alt blev 54 celletyper, for hvilke både CTCF ChIP-seq-data og RNA-seq-data er offentligt tilgængelige, udvalgt (Additional file 6: Table S5) til undersøgelse af sammenhængen mellem CTCF-binding og genekspression for hvert CTCF-genpar på det samme kromosom. For at opnå CTCF-bindingsniveauet blev der genereret en læsningstællingsmatrix ved hjælp af læsninger pr. kilobase pr. million (RPKM) på union CTCF-bindingssteder fra ChIP-seq-data. Læsningstællingsmatrixen blev skaleret med kvadratroden af RPKM efterfulgt af kvantilnormalisering. Genekspressionsniveauet blev målt for hvert gen ved hjælp af kvadratroden af transskriptioner pr. million (TPM) fra RNA-seq-data. For hvert CTCF-genpar kvantificerede vi forbindelsen mellem CTCF-stedet og genet på tværs af alle 54 celletyper ved hjælp af korrelationskoefficienten R mellem det normaliserede CTCF-bindingsniveau og genekspressionen (fig. 3a). CTCF-genpar blev anset for at være “stærkt korreleret” med R2 større end 0,25, f.eks, korrelationskoefficient større end 0,5 eller mindre end – 0,5, og de stærkt korrelerede CTCF-genpar bidrager til 1,3 % af alle CTCF-genpar (Additional file 1: Fig. S8a).
Identifikation af konstitutive CTCF-bundne kromatindomæner
For hvert CTCF-bindingssted definerede vi dets tilknyttede kromatindomæne som den genomiske region, der (1) omfatter dette specifikke CTCF-bindingssted, (2) er afgrænset af et par konstitutive CTCF-bindingssteder med motiver af modsatte orienteringer, og (3) optager mindst 100 kb og højst 1 MB-region på hver side af CTCF-bindingsstedet. Figur 3b indeholder en skematisk oversigt over, hvordan konstitutive CTCF-bundne kromatin-domæner blev defineret.
Detektion af DNA-methyleringsændringer omkring CTCF-bindingssteder
DNA-methyleringsændringer blev detekteret inden for en 300-bp-region centreret ved hvert CTCF-bindingssted. Regioner med mindst 3 CpG’er dækket af mindst 5 læsninger (≥ 5×) i både kræftcellelinjer og tilsvarende normale væv blev bevaret. En 300-bp-region blev detekteret som differentielt methyleret, hvis de gennemsnitlige differentielle methyleringsniveauer for alle CpG’er (≥ 5×) inden for denne region var større end 20 % .
Detektion af mutationshastighed og differentiel motivscore
For hvert CTCF-bindingssted blev det rå mutationstallet beregnet som forekomsten af mutationshændelser i alle prøver/patienter ved hvert enkelt basepar inden for en 400-bp-region centreret på CTCF-bindingsstedet. Mutationsraten for en gruppe af CTCF-bindingssteder blev beregnet som det gennemsnitlige antal mutationer over antallet af CTCF-bindingssteder for hvert basepar inden for 400-bp-regionen.
Motivscore blev målt ved at score CTCF-positionsvægtmatrixen (Jaspar , Matrix ID: MA0139.1) til en 19-bp DNA-sekvens centreret ved CTCF-motivet eller CTCF-bindingsstedet ved hjælp af log sandsynlighedsforhold (med baggrundsnukleotidfrekvens som for A,C,G,T). Den differentielle motivscore blev beregnet ved at sammenligne motivscorerne for reference- og de muterede sekvenser.
DNA-sekvensmotivanalyse
DNA-sekvensmotivberigelsesanalyse blev udført ved hjælp af MDSeqPos (version 1.0.0.0) på Cistrome med standardparametre (-cisrome -Homo Sapien eller Mus musculus). De novo-motivanalyser blev udført ved hjælp af HOMER (version 4.10) med findmotifs.pl-modulet og MEME (version 5.1.1) med følgende parametre: meme -dna -mod zoops -maxw 20 -evt -0.01.
Identifikation af CTCF intra-domæne differentielt interagerede regioner
For et givet sæt af CTCF-bindingssteder blev kromatininteraktionsændringerne mellem et CTCF-sted og hver af dets intra-domæne ikke-overlappede bins, målt fra normaliserede Hi-C-kontaktkort i kræftceller i forhold til matchede normale celler, indsamlet for hvert af CTCF-bindingsstederne (Additional file 1: Fig. S14b). Regioner med nedsatte interaktioner (log2 FC < -1, gennemsnitlig log2-interaktion > 0) med kræftspecifikke tabte CTCF-bindingssteder, og regioner med øgede interaktioner (log2 FC > 1, gennemsnitlig log2-interaktion > 0) med kræftspecifikke vundne CTCF-bindingssteder blev brugt til downstream transkriptionsfaktor (TF)-berigelsesanalyse.
Transkriptionsfaktorberigelsesanalyse
En revideret version af BART-algoritmen blev anvendt til TF-berigelsesanalyse. Kort fortalt blev en samling af union DNase I hypersensitive steder (UDHS) tidligere kurateret som et repertoire af alle kandidat-cis-regulerende elementer i det menneskelige genom, og 7032 ChIP-seq-datasæt blev indsamlet for 883 TF’er , hvor hver TF havde et eller flere ChIP-seq-datasæt fra flere celletyper eller betingelser. Der blev genereret en binær profil for hver TF på UDHS, der angiver, om TF’en har mindst én top fra et af sine ChIP-seq-datasæt placeret inden for hver af UDHS’erne. Der blev anvendt en analyse af bindingsberigelse for hver TF ved at sammenligne TF-bindingen på en delmængde af UDHS, der overlapper de udvalgte genomiske regioner, med TF-bindingen på UDHS. p-værdien blev opnået ved hjælp af Fisher’s eksakte test.