Experimentellt tillvägagångssätt
Patientens xenografting och cellodling
De humana T-ALL-cellinjerna inkluderar CUTLL1 (gåva från Adolfo Ferrando, Columbia University) och JURKAT (American Type Culture Collection (ATCC), Manassas, VA, #CCL-119). Cellerna odlades i RPMI1640-medium med l-glutamin och 25 mM HEPES (Corning) kompletterat med 10 % värmeinaktiverat fetalt bovint serum (Sigma-Aldrich), 10 U/mL penicillin-streptomycin (Gibco) och 1× glutaMAX (Gibco) i en luftfuktad inkubator vid 37 °C och 5 % CO2. Cellerna testas regelbundet för förekomst av mykoplasma med hjälp av Lonza Walkersville MycoAlert Mycoplasma Detection Kit (sista testet i januari 2020). Cellinjerna hålls i kultur i högst 20 passager och autentiseras med hjälp av profilering av korttandemrepeater (JURKAT) eller med hjälp av PCR för att detektera TCRb-NOTCH1-translokationen (TCRBJ2S4CUTLL1F:5′-GGACCCGGCTCTCAGAGTGCT-3′, NOTCH1CUTTL1R:5′-TCCCGCCCTCTCCAAAATAAGG-3′). Den sista cellverifieringen utfördes i februari 2020. Mänskliga CD4+ T-celler köptes från AllCells. Primära mänskliga prover samlades in med informerat samtycke och analyserades under överinseende av institutionella granskningsnämnden vid Padovas universitet, Associazone Italiana di Ematologia e Oncologia Pediatrica och Berlin-Frankfurt-Münster (AIEOP-BFM) ALL 2000/2006 pediatric clinical trials. Informerat samtycke till att använda överblivet material för forskningsändamål inhämtades från alla patienter när de gick in i studien i enlighet med Helsingforsdeklarationen.
Antikroppar och reagenser
Western blots utfördes med följande antikroppar: Actin och CTCF från Millipore Sigma (klon C4; 07-729) och kluven NOTCH1 (Val1744) från Cell Signaling Technology (4147). ChIP-seq utfördes med följande antikroppar: CTCF från Millipore Sigma (07-729), H3K27Ac (8173S) och H3K27me3 (9733S) från Cell Signaling Technology och H3K4me1 (07-473) från Millipore.
In situ Hi-C
In situ Hi-C utfördes på CD4+ T-celler, Jurkat, CUTLL1 och patientens xenografts enligt tidigare beskrivningar . I korthet korslänkades cellerna med 1 % formaldehyd i 10 minuter i rumstemperatur. Per Hi-C-reaktion lyserades 5 miljoner celler och kärnorna permeabiliserades. DNA smältes med MboI från New England Biolabs (R0147M). De uppdelade fragmenten märktes med biotinylerat d-ATP från Jena Bioscience (NU-835-BIO14-S) och ligerades. Efter RNase-behandling och Proteinase K-behandling för att upphäva tvärbindningar sonikades kärnorna med hjälp av en Covaris E220 för att producera en genomsnittlig fragmentlängd på 400 bp. Streptavidinpärlor från Thermo Fisher Scientific (65001) användes för att dra ner biotinmärkta fragment. Efter rening och isolering av DNA framställdes slutliga bibliotek med hjälp av NEBNext® Ultra™ II DNA Library Prep Kit for Illumina® och sekvenserades med hjälp av parvis sekvensering vid en avläsningslängd på 150 bp på en Illumina HiSeq 2500 för att producera i genomsnitt 400 miljoner avläsningar per prov.
ChIP-seq-profilering
CD4+ T-celler, Jurkat, CUTLL1 och patientens xenografts korslänkades med 1 % formaldehyd och 1 % fetalt bovint serum i PBS i 10 minuter i rumstemperatur. Reaktionen släcktes med 0,2 M glycin i rumstemperatur i 5 minuter. Cellerna tvättades sedan med PBS och pelleterades.
För CTCF ChIPs utfördes immunutfällning enligt ett protokoll som beskrivits tidigare . En pellet innehållande 50 miljoner celler lyserades med 5 ml lysisbuffert (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 minuter vid 4 °C. Kärnorna pelleterades vid 1350×g i 7 minuter och resuspenderades i 10 mM Tris pH 8, 1 mM EDTA och 0,1 % SDS. Kromatinet klipptes med ett Covaris E220-system till en genomsnittlig fragmentlängd på 400 bp och snurrades vid 15 000 rpm i 10 minuter för att avlägsna olösligt kromatin och skräp. Supernatanten inkuberades med 20 μL Dynabeads Protein G i 30 minuter innan pärlorna kasserades. En procent av den totala volymen sparades som input och resten inkuberades med anti-CTCF-antikropp över natten. Totalt tillsattes 100 μL Dynabeads Protein G under 2 h. De bundna fragmenten tvättades två gånger med 1 mL lågsaltsbuffert (20 mM Tris-HCl pH 8,0, 150 mM NaCl, 2 mM EDTA, 1 % w/v Triton X-100 och 0.1 % w/v SDS), en gång med högsaltsbuffert (20 mM Tris-HCl pH 8,0, 500 mM NaCl, 2 mM EDTA, 1 % w/v Triton X-100 och 0,1 % w/v SDS), en gång med litiumkloridbuffert (10 mM Tris-HCl pH 8.0, 250 mM LiCl, 1 mM EDTA, 1 % w/v NP-40 och 1 % w/v deoxycholsyra) och två gånger med TE (10 mM Tris pH 8, 1 mM EDTA).
För histon ChIP:s lyserades cellerna i 375 μL av en inkubationsbuffert för kärnor (15 mM Tris pH 7.5, 60 mM KCl, 150 mM NaCl, 15 mM MgCl2, 1 mM CaCl2, 250 mM sackaros, 0,3 % NP-40, 1 mM NaV, 1 mM NaF och 1 EDTA-fri proteashämmartablett (Roche)/10 mL i H2O) i 10 minuter på is. Kärnorna tvättades en gång med digestbuffert (10 mM NaCl, 10 mM Tris pH 7,5, 3 mM MgCl2, 1 mM CaCl2, 1 mM NaV, 1 mM NaF och 1 EDTA-fri proteashämmartablett (Roche)/10 mL i H2O) och resuspenderades i 57-μL digestbuffert som innehåller 4,5 enheter MNas (USB) i 1 h vid 37 °C. MNaseaktiviteten släcktes i 10 minuter på is efter tillsats av EDTA till en slutkoncentration på 20 mM. Kärnorna pelleterades och resuspenderades i 300-μL nukleolysbuffert (50 mM Tris-HCl pH 8,0, 10 mM EDTA pH 8,0, 1 % SDS, 1 mM NaV, 1 mM NaF och 1 EDTA-fri proteashämmare i tablettform (Roche)/10 mL i H2O) före sonikation med en Bioruptor Pico (Diagenode) i 5 min (30 s på, 30 s av). Lysatet centrifugerades vid maxhastighet i 5 minuter för att avlägsna skräp. Nio volymer IP-spädningsbuffert (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 och 1 EDTA-fri proteashämmare i tablettform (Roche)/10 mL i H2O) tillsattes till supernatanten. Totalt tillsattes 50 μL Dynabeads Protein G och provet inkuberades vid 4 °C i 30 minuter med rotation. En procent av provet behölls som input och det återstående provet delades upp i tre rör. Totalt 50 μL Dynabeads Protein G konjugerat med 15 μL av lämplig antikropp tillsattes till varje rör före inkubation över natten vid 4 °C, roterande. Pärlbundna komplex tvättades i 5 minuter vardera i 1 ml lågsaltbuffert, högsaltbuffert, LiCl-buffert och två gånger med TE.
För att eluera pärlbundna komplex resuspenderades pärlorna i 50 μL elutionsbuffert (100 mM NaHCO3, 1 % w/v SDS) och inkuberades vid 65 °C i 15 minuter, med skakningar vid 1 000 varv per minut på en termomixer (Thermo Scientific). Elueringen upprepades en andra gång och sedan tillsattes 100 μL RNasbuffert (12 μL 5 M NaCl, 0,2 μL 30 mg/mL RNas och 88 μL TE) till varje ChIP- och inmatningsprov. Proverna inkuberades vid 37 °C i 20 minuter, följt av tillsats av 100 μL proteinas K-buffert (2,5 μL 20 mg/mL proteinas K, 5 μL 20 % SDS och 92,5 μL TE) över natten vid 65 °C. En lika stor volym fenol:kloroformlösning tillsattes och blandades noggrant. Blandningen överfördes till MaXtract High Density-rör (Qiagen) och centrifugerades i 8 minuter vid 15 000 rpm. Den övre fasen överfördes till nya rör och blandades med 1,5 μL 20 mg/mL glykogen, 30 μL 3M natriumacetat och 800 μL etanol. Proverna inkuberades vid – 80 °C tills de frystes och centrifugerades sedan vid 15 000 rpm i 30 minuter vid 4 °C. Supernatanten avlägsnades och pellets tvättades i 800 μL 70 % iskall etanol och centrifugerades i 10 minuter vid 4 °C och 15 000 rpm. Efter noggrant avlägsnande av etanolen lufttorkades pellets och resuspenderades i 30 μl 10 mM Tris vid pH 8.
IP och ingående DNA kvantifierades sedan med hjälp av en Qubit 3.0-fluorometer. Bibliotek framställdes med hjälp av KAPA HyperPrep Kit (KK8505) och sekvenserades med en Illumina NextSeq 500 till ett genomsnittligt djup på 28 miljoner reads per prov.
RNA-seq-profilering
RNA isolerades från 3 miljoner celler per prov med hjälp av Bio-Rad Aurum™ Total RNA Mini Kit och kvantifierades med Agilent RNA 6000 Nano Kit med Agilent Bioanalyzer. Bibliotek framställdes genom rRNA-depletion med hjälp av Illumina TruSeq® Stranded mRNA Library Prep Kit för en låg koncentration av startprovet och sekvenserades genom single end-sekvensering på en Illumina NextSeq 500 till ett genomsnittligt djup på 18 miljoner reads per prov.
DNA-metyleringsprofilering
Genomiskt DNA isolerades med hjälp av AllPrep DNA/RNA Micro Kit (Qiagen). För att bedöma DNA-metyleringsstatus för hela genomet utförde vi mRRBS . Efter fluorometrisk kvantifiering med hjälp av ett Qubit 3.0-instrument smälte vi genomiskt DNA med restriktionsenzymet MspI (New England Biolabs) och storlekssorterade för fragment som var ungefär 100-250 baspar långa med hjälp av SPRI-pärlor (solid phase reversible immobilization) (MagBio Genomics). Det resulterande DNA:t genomgick bisulfitkonvertering med hjälp av EZ DNA Methylation-Lightning Kit (Zymo Research). Vi skapade bibliotek från bisulfitkonverterat enkelsträngat DNA med hjälp av Pico Methyl-Seq Library Prep Kit (Zymo Research), som sedan sammanfördes för sekvensering på ett Illumina NextSeq 500-instrument med hjälp av NextSeq 500/550 V2 High Output-reagenskit (1 × 75 cykler) till ett minsta läsdjup på 50 miljoner läsningar per prov.
Helgenomsekvensering
Tre miljoner celler från cellinjer eller patientprover pelleterades och resuspenderades i 1 mL Cell Lysis Solution (Qiagen) blandat med 500 μg RNas A. Lysireaktionen genomfördes vid 37 °C i 15 minuter. Totalt tillsattes 333 μL proteinutfällningslösning (Qiagen) till varje prov som sedan vortexades och därefter centrifugerades vid 2000×g i 10 min. Supernatanten blandades med 1 mL isopropanol tills DNA-strängar fälldes från lösningen. Efter att ha kastat supernatanten tvättades DNA-pelleten med 1 mL 70 % etanol och centrifugerades vid 2000×g i 1 min. Etanolen hälldes sedan ut och pelleten lufttorkades i 15 minuter innan den resuspensionerades i 50-100 μL DNA-hydreringslösning (Qiagen). DNA sekvenserades med parvis Illumina-sekvensering med 30× täckning.
Immunutfällning
Totalt 100 miljoner celler för varje immunutfällningsreaktion pelleterades och inkuberades i buffert A (10 mM HEPES pH 8,0, 1,5 mM MgCl2, 10 mM KCl, 0,5 mM DTT) i 10 min på is. Cellerna lyserades sedan med 12 slag med en 7 ml vävnadskvarn med lös pistill (Wheaton, 357542) och centrifugerades vid 2000 rpm i 7 min. Kärnpellets resuspenderades i 5 volymer TENT-buffert (50 mM Tris pH 7,5, 5 mM EDTA, 150 mM NaCl, 1 % Triton X-100, 5 mM MgCl2) och behandlades med bensonas i 30 min före 5 passager genom en 25 g × 5/8 tums spruta. Den olösliga fraktionen avlägsnades efter centrifugering vid 2000 rpm i 7 minuter och inkuberades över natten med Dynabeads Protein G hybridiserad med antikropp. Sammanlagt 2 miljoner celler avlägsnades för input. Lysat av pärlor och kärnor tvättades 6 gånger med TENT-buffert och eluerades sedan i 0,1 M glycin pH 2,5 med 100 mM Tris pH 8,0 före. NuPAGE LDS-provbuffert tillsattes till eluaten och inmatningarna, som sedan inkuberades vid 70 °C i 15 minuter före analys med Western blot.
Offentlig datainsamling
Offentliga CTCF ChIP-seq-data samlades in från Cistrome Data Browser (för toppfiler) och NCBI GEO (för fastq-filer, Additional file 2: Table S1). ChIP-seq-data om histonmodifiering samlades in från NCBI GEO och ENCODE (för bam-filer). Offentliga RNA-seq-data i flera celltyper samlades in från ENCODE (för fastq-filer). Data om DNA-metyleringsprofilering samlades in från ENCODE (för bed bedMethyl-filer) och NCBI GEO. Hi-C-data samlades in från NCBI GEO och ENCODE (för fastq-filer). ATAC-seq-data samlades in från NCBI GEO (för fastq-filer). Helgenomsekvenseringsdata för BRCA-, COAD-, LUAD- och PRAD-prover samlades in från International Cancer Genome Consortium (ICGC) Data Portal . Detaljerad information inklusive accessions-id för alla offentliga dataset som samlats in i detta arbete finns i Additional file 6: Table S5.
Databehandling
Analys av ChIP-seq-data
Sekvensanpassning för ChIP-seq-data i fastq-filer utfördes med hjälp av samma standardiserade analyspipeline som används i Cistrome DB , för att uppnå konsistens och reproducerbarhet. Alla sekvensdata genomisk anpassning utfördes med hjälp av Chilin-pipeline med standardparametrar ($ chilin simple -p narrow -s hg38 –threads 8 -t IN.fq -i PRENAME -o OUTDIR). I korthet anpassades sekvensläsningarna till det mänskliga referensgenomet (GRCH38/hg38) med hjälp av 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-filerna konverterades sedan till bam-filer med hjälp av samtools ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). För CTCF ChIP-seq-datasetterna användes MACS2 för att kalla toppar under FDR-tröskelvärdet 0,01 ($ macs2 callpeak –SPMR -B -q 0,01 –keep-dup 1 -g hs -t PRENAME.bam -n PRENAME –outidr OUTDIR). Toppar med en viktad anrikning på minst 4 bibehölls. Bigwiggle-filer genererades med hjälp av BEDTools och UCSC-verktyg ($ 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). Slutligen inkluderades endast de CTCF ChIP-seq-prover som har minst 2 000 toppar i den integrativa analysen nedströms.
ATAC-seq-dataanalys
Trim Galore användes för att trimma de obearbetade sekvenseringsläsningarna ($ trim_galore –nextera –phred33 –fastqc –paired R1.fq R2.fq -o OUTDIR). Läsningarna anpassades till det mänskliga referensgenomet (GRCH38/hg38) med hjälp av Bowtie2 ($ bowtie2 -p 10 -X 2000 -x INDEX -1 R1.fq -2 R2.fq -S PRENAME.sam). Sam-filerna konverterades sedan till bam-filer med hjälp av samtools ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). Bedtools användes för att konvertera bam-filer till bed-format ($ bamToBed -i PRENAME.bam -bedpe > PRENAME_PE.bed). Läsningar som kartlades till mitokondrie-DNA kasserades från nedströmsanalysen.
RNA-seq-dataanalys
RNA-seq-datasatser bearbetades med hjälp av Salmon ($ salmon quant –gcBias -i INDEX -l A -p 8 {-1 R1.fq -2 R2.fq| -r IN.fq} -o OUTDIR). Transkriptomindex byggdes på det mänskliga referensgenomet (GRCH38/hg38). Uppskattningar av abundansen på transkriptnivå sammanfattades till gennivå med hjälp av paketet ”tximport” för differentiell uttrycksanalys. DESeq2 användes för att identifiera differentiellt uttryckta gener, och olika tröskelvärden som användes i olika analyser angavs på motsvarande sätt i manuskriptet.
Hi-C dataanalys
Hi-C data bearbetades med hjälp av HiC-Pro ($ HiC-Pro -i INDIR -o OUTDIR -c CONFIG -p). Kontaktkartor genererades med en upplösning på 5 kb. Rå matrisdata normaliserades med hjälp av den metod som beskrivs i Normalization of Chromatin Interactions.
DNA-metyleringsdataanalys
DNA-metyleringsdata (för T-ALL-cellinjer och T-ALL-patienter) demultiplexerades med bcl2fastq följt av trimmning av 10 baspar från 5′-ändan för att ta bort primer- och adaptersekvenser med TrimGalore . Sekvensanpassning till referensgenomet GRCh38/hg38 och metyleringskallelser utfördes med Bismark ($ bismark –multicore 8 –bowtie2 -q -N 1 INDEX INFILE.fq). Täckningsfiler (counts) för cytosiner i CpG-kontext genererades med Bismark ($ bismark_methylation_extractor –multicore 8 –comprehensive –bedGraph INFILE_bismark_bt2.bam).
Analys av sekvenseringsdata för hela genomet
Mutationer identifierades för två T-ALL-cellinjer (Jurkat och CUTLL1) och två T-ALL-patientprover från sekvenseringsdata för hela genomet. Vi anpassade Illumina short-read-sekvenserna till det mänskliga referensgenomet (GRCH38/hg38) med hjälp av BWA mem. Vi använde SAMBlaster för att identifiera de disharmoniska paren, dela upp läsningar och flagga de förmodade PCR-duplikaten. Vi använde SAMBAMBA för att konvertera SAM-justeringen till BAM-format, och samtools användes för att sortera de justerade för att skapa en BAM-fil som motsvarar varje prov.
Vi använde VarDict för att identifiera de varianter som överlappade unionens CTCF-bindningsställen. Vi använde alla standardparametrar utom ”-f 0,1” som användes för att identifiera varianter som stöddes av mer än 10 % av läsningarna på den platsen. Vi annoterade varianterna med hjälp av Variant Effect Predictor (VEP) och använde anpassade skript för att identifiera de varianter som påverkar TF-bindningen.
Vi använde återigen VarDict för att identifiera varianterna i CTCF- och NOTCH1-generna för de fyra proven. Vi använde alla standardparametrar utom ”-f 0,1” som användes för att identifiera varianter som stöddes av mer än 10 % av läsningarna på den platsen. Vi annoterade varianterna med hjälp av Variant Effect Predictor (VEP) , och filtrerade den sedan för att identifiera de mutationer som antingen a) inte förekom i mer än 1 % av en normal mänsklig population, eller b) hade en CADD-poäng för deleteriousness > 20, eller c) fanns i COSMIC-databasen.
Integrativ modellering och statistisk analys
Identifiering av CTCF-bindningsrepertoaren i det mänskliga genomet
För CTCF ChIP-seq samlade vi in totalt 793 dataset, inklusive 787 offentliga dataset och 6 dataset som vi genererade (Additional file 2: Table S1). Totalt användes 771 CTCF ChIP-seq-dataset med toppar över 2000 i den här studien. Varje dataset kan ge MACS2-identifierade CTCF-toppar i intervallet mellan 2050 och 198 021, med en median på 46 451 och totalt 36 873 077 toppar (Additional file 1: Fig. S1a). Fördelningen av intervalllängderna mellan intilliggande toppar av CTCF-toppar för alla 36 873 077 toppar från de 771 datasetten har en böjningspunkt vid ~ 150 bp (Additional file 1: Fig. S1c) som anger gränsen mellan samma bindningsställe och olika bindningsställen . Därför använde vi 150 bp som gränsvärde för att slå samman CTCF-toppar. I praktiken förlängde vi ± 75 bps från varje topp för att generera en 150 bp-region centrerad vid toppen för att representera varje topp och slog samman alla överlappande toppregioner för att generera en unionsuppsättning av CTCF-bindningsställen, som innehåller 688 429 icke överlappande ställen. Varje bindningsplats tilldelades en CTCF-beläggningspoäng, definierad som antalet ChIP-seq-datasatser som uppvisar en topp inom platsen. Följaktligen definierade vi ockupationsfrekvensen som förhållandet mellan ockupationsvärdet och det totala antalet CTCF ChIP-seq-datasetterna. För att ytterligare säkerställa robustheten hos de identifierade CTCF-bindningsplatserna valde vi ut 285 467 platser med hög konfidens med ockupationspoäng ≥ 3 för nedströmsanalyser. CTCF-motiv inom de fackliga bindningsplatserna söktes med FIMO med Jaspar-matris (ID: MA0139.1), med en p-värdeströskel på 1e-4. Ett motiv med det minsta p-värdet behölls för varje CTCF-bindningsplats.
Identifiering av konstitutiva CTCF-bindningsplatser
Fördelningen av ockupationspoängen för alla 285 467 CTCF-bindningsplatser (Additional file 1: Fig. S1d, blå kurva) visar att majoriteten av CTCF-bindningsplatserna förekommer i endast ett fåtal dataset, och att antalet bindningsplatser minskar med ökande ockupationspoäng när ockupationspoängen är liten. Det finns dock CTCF-bindningsställen som är mycket bevarade i nästan alla dataset (t.ex. bindningsställen med en beläggningspoäng på mer än 600). Vi använder en power law-funktion för att anpassa distributionskurvan (blå) som visas i Additional file 1: Fig. S1d för att bestämma gränsvärdet för konstitutiva CTCF-platser. Vi betecknar Oi som antalet observerade CTCF-bindningsställen med en beläggningspoäng som är lika med i, och Ei som antalet förväntade CTCF-ställen med en beläggningspoäng som är lika med i. Anpassningen av kraftlagen till data Oi kan beskrivas som (Additional file 1: Fig. S1d, grönt):
Vi definierar gränsvärdet A för konstitutiva CTCF-bindningsställen som:
Med andra ord bör de totala observerade CTCF-platserna med en beläggningspoäng som är större än A vara 6 gånger fler än förväntat. Vi fastställde sedan A = 615 och använde en gräns för ockupationsfrekvens på 80 % för att definiera 22 097 konstitutiva CTCF-bindningsställen, vilket motsvarar ockupationspoängen ≥ 616 i alla 771 CTCF ChIP-seq-dataset.
Identifiering av cancerspecifika vunna/förlorade CTCF-bindningsställen
Vi använde följande 2 kriterier för att identifiera cancerspecifika förlorade CTCF-bindningsställen: (1) CTCF-bindningsplatsen ska ha en lägre beläggningsfrekvens för dataset för den cancertypen jämfört med beläggningsfrekvensen för alla dataset och (2) CTCF-bindningsnivån (kvantifierad som normaliserade ChIP-seq-läsningssiffror) vid platsen är lägre i cancerdatasetterna än i andra dataset. För vunna CTCF-platser använde vi den omvända uppsättningen kriterier. Kort sagt, för varje CTCF-bindningsplats i varje cancertyp beräknades ockupationsvärdet i cancerdatasetterna tillsammans med dess ockupationsvärde i alla 771 dataset. CTCF-bindningsnivåerna erhölls från en normaliserad matris för antalet läsningar där ChIP-seq-läsningssiffrorna (RPKM) först beräknades för fackliga CTCF-bindningsställen i alla dataset och sedan följdes av en kvantilnormalisering. Vi använde opartade tvåsidiga Student’s t-test för att kvantifiera skillnaden i bindningsnivåer mellan olika grupper av dataset, och p-värdet justerades sedan med hjälp av Benjamini-Hochberg-proceduren . Dessutom jämfördes bindningsbeläggningspoäng och bindningsnivåer mellan dataset för cancer och dataset från matchade normala vävnader eller celltyper, för att ta hänsyn till den potentiella förväxlingsfaktorn vävnadsspecificitet snarare än cancerspecificitet. Detaljerade kriterier för att identifiera cancerspecifika CTCF-bindningsställen beskrivs nedan:
-
Cancerspecifika förlorade CTCF-bindningsställen: (1) beläggningsfrekvens ≤ 0,2 i cancerdatabaser; (2) beläggningsfrekvens ≥ 0,7 i 771-databaser; (3) beläggningsfrekvens ≥ 0,2 i cancerdatabaser; (4) beläggningsfrekvens ≥ 0,5 i 771-databaser.5 (med beläggningsgrad ≥ 2) i matchade dataset för normal vävnad; (4) CTCF-nivåerna är lägre i cancer jämfört med alla andra dataset (statistisk poäng < 0), (5) CTCF-nivåerna är lägre i cancer jämfört med matchade dataset för normal vävnad (statistisk poäng < 0), (6) genomsnittliga CTCF-bindningssignaler (RPKM) < 5 i cancerdataset.
-
Cancerspecifika vunna CTCF-bindningsställen: (1) beläggningsfrekvens ≥ 0,5 (med beläggningspoäng ≥ 2) i cancerdataset, (2) beläggningsfrekvens ≤ 0,2 i 771 dataset, (3) beläggningspoäng = 0 i matchade dataset för normal vävnad, (4) CTCF-nivåer är signifikant högre i cancer jämfört med alla andra dataset (FDR ≤ 0.01), (5) CTCF-bindningsnivåerna är signifikant högre i cancer jämfört med matchade dataset för normal vävnad (FDR ≤ 0,01), (6) genomsnittliga CTCF-bindningssignaler (RPKM) > 2 i dataset för cancer.
De specifika vunna och förlorade CTCF-bindningsställena för varje cancertyp visas i Additional file 4: Table S3.
Kvantifiering av differentiell kromatintillgänglighet
Vi använde de bearbetade data från Ref. som inkluderar en matris av normaliserade ATAC-seq-insättningssiffror inom TCGA-pan-cancer-toppsetet för att bedöma den differentiella kromatintillgängligheten runt CTCF-bindningsställena. För varje cancertyp bland BRCA, CRC, LUAD och PRAD användes de pan-cancer ATAC-seq-toppar som överlappar med identifierade cancertypspecifika förlorade eller ökade CTCF-bindningsställen för nedströmsanalyser. ATAC-seq-differentialpoängen för varje topp kvantifierades som vikförändringen av genomsnittet av de normaliserade ATAC-seq-insättningssiffrorna från patientprover i motsvarande cancertyp jämfört med patienter i andra cancertyper, och ATAC-seq-differentialpoängen tilldelades sedan toppens överlappande CTCF-bindningsställe.
För att få konsekvens tillämpade vi samma tillvägagångssätt som används för TCGA ATAC-seq-data för att analysera de insamlade ATAC-seq-data från T-ALL-cellinjen Jurkat och normala CD4+ T-celler. En datamatris genererades med hjälp av ATAC-seqs råa avläsningssiffror på fackliga CTCF-bindningsställen för alla dataset för Jurkat och T-celler. Kvantelnormalisering tillämpades på den log2-skalade matrisen (pseudoantal = 5). ATAC-seq-differentialpoängen mättes som vikförändringen av de genomsnittliga normaliserade ATAC-seq-räkningarna mellan dataset för Jurkat respektive CD4+ T-celler vid varje CTCF-bindningsplats.
Normalisering av kromatininteraktioner
Givet en Hi-C-kontaktkarta A = {aij} återspeglar poängen aij kartlagda läsningar mellan två genomiska regioner i och j. Antag att binstorleken är 5 kb, kommer regionerna i och j att ha ett genomiskt avstånd på ∣i – j ∣ × 5kb. Eftersom kontaktsannolikheten mellan två bin minskar med ökande genomiskt avstånd , normaliserade vi kontaktkartan på följande sätt: För varje givet genomiskt avstånd dk = k × 5kb kvantifierar vi en normaliseringsfaktor \( { {\overline{S}}_{d_k} \) som de genomsnittliga interaktionerna mellan alla bin-par med samma genomiska avstånd dk i en och samma kromosom, t.ex, \( {\overline{S}}_{d_k}=\left({\sum}_{j-i=k}{a}_{ij}\right)/n \), där n är det totala antalet bin-par med avstånd dk. Interaktionspoängen aij mellan två bin med avstånd dk normaliserades sedan med \( {\overline{S}}_{d_k} \) som \( {a}_{ij}^{\prime }={a}_{ij}/{\overline{S}}_{d_k} \). Med detta tillvägagångssätt normaliserade vi matrisen A till \( A^{\prime }=\left\{{a}_{ij}^{\prime}\right\} \) inom varje kromosom.
Detektering av differentiella kromatininteraktioner
Vi betecknade de normaliserade Hi-C-kontaktkartorna i cancerdataset och det normala datasetetet som C = {cij} respektive N = {nij}. För en given CTCF-bindningsplats x (med koordinat xc) och ett fördefinierat genomiskt avstånd L samlas kromatininteraktionerna mellan x och dess närliggande icke överlappande 5-kb-bins med genomiskt avstånd upp till L in från C respektive N. Interaktionspoäng mellan x och dess närliggande 5-kb-bins i C samlas in som IC = {cij} , medan antingen i eller j är lika med ⌊xc/5kb⌋ och 0 < (j – i) × 5kb ≤ L. På samma sätt samlades interaktionspoängen mellan x och dess närliggande 5-kb-bins i N som IN = {nij}. Ett parat tvåsidigt Student’s t-test tillämpades sedan på IC och IN för att kvantifiera den differentiella interaktionen mellan cancerceller och normala celler som omger CTCF-bindningsstället x.
Samband mellan CTCF-bindning och genuttryck
Totalt valdes 54 celltyper för vilka både CTCF ChIP-seq-data och RNA-seq-data är offentligt tillgängliga (Additional file 6: Table S5) för att undersöka sambandet mellan CTCF-bindning och genuttryck för varje CTCF-genpar på samma kromosom. För att få fram CTCF-bindningsnivån genererades en matris för antalet läsningar med hjälp av läsningar per kilobas per miljon (RPKM) på fackliga CTCF-bindningsställen från ChIP-seq-data. Läsräkningsmatrisen skalades med kvadratroten av RPKM följt av kvantilnormalisering. Genuttrycksnivån mättes för varje gen med hjälp av kvadratroten av transkriptioner per miljon (TPM) från RNA-seq-data. För varje CTCF-genpar kvantifierade vi sambandet mellan CTCF-platsen och genen över alla 54 celltyper med hjälp av korrelationskoefficienten R mellan den normaliserade CTCF-bindningsnivån och genuttrycket (fig. 3a). CTCF-genpar bedömdes som ”starkt korrelerade” med R2 större än 0,25, t.ex, korrelationskoefficient större än 0,5 eller mindre än – 0,5, och de högt korrelerade CTCF-genparen bidrar till 1,3 % av alla CTCF-genpar (Additional file 1: Fig. S8a).
Identifiering av konstitutiva CTCF-bundna kromatindomäner
För varje CTCF-bindningsställe definierade vi dess associerade kromatindomän som den genomiska region som (1) innefattar detta specifika CTCF-bindningsställe, (2) begränsas av ett par konstitutiva CTCF-bindningsställen med motiv med motsatt inriktning, och (3) upptar minst 100 kb och högst 1 MB-region på vardera sidan av CTCF-bindningsstället. Figur 3b innehåller en schematisk bild av hur konstitutiva CTCF-bundna kromatindomäner definierades.
Detektering av DNA-metyleringsförändringar som omger CTCF-bindningsställen
DNA-metyleringsförändringar detekterades inom en 300-bp-region centrerad vid varje CTCF-bindningsställe. Regioner med minst 3 CpGs som täcktes av minst 5 läsningar (≥ 5×) i både cancercelllinjer och motsvarande normala vävnader behölls. En 300-bp-region upptäcktes som differentiellt metylerad om de genomsnittliga differentiella metyleringsnivåerna för alla CpGs (≥ 5×) inom denna region var större än 20 % .
Detektering av mutationsfrekvens och differentiell motivpoäng
För varje CTCF-bindningsställe beräknades det obearbetade mutationsantalet som förekomsten av mutationshändelser i alla prover/patienter vid varje enskilt baspar inom en 400-bp-region centrerad vid CTCF-bindningsstället. Mutationsfrekvensen för en grupp av CTCF-bindningsställen beräknades som det genomsnittliga antalet mutationer över antalet CTCF-bindningsställen för varje baspar inom 400-bp-regionen.
Motivpoäng mättes genom att poängsätta CTCF-positionsviktsmatrisen (Jaspar , Matrix ID: MA0139.1) för en 19-bp-dna-sekvens centrerad vid CTCF-motivet eller CTCF-bindningsstället med hjälp av log sannolikhetskvoter (med bakgrundsnukleotidfrekvensen som för A,C,G,T). Den differentiella motivpoängen beräknades genom att jämföra motivpoängen för referenssekvenserna och de muterade sekvenserna.
DNA-sekvensmotivanalys
DNA-sekvensmotivberikningsanalysen utfördes med hjälp av MDSeqPos (version 1.0.0) på Cistrome med standardparametrar (-cisrome -Homo Sapien eller Mus musculus). De novo-motivanalyser utfördes med hjälp av HOMER (version 4.10) med modulen findmotifs.pl och MEME (version 5.1.1) med följande parametrar: meme -dna -mod zoops -maxw 20 -evt -0,01.
Identifiering av CTCF:s intradomäner med differentiell interaktion
För en given uppsättning CTCF-bindningsställen samlades kromatininteraktionsförändringarna mellan en CTCF-plats och var och en av dess intradomäner utan överlappning, uppmätta från normaliserade Hi-C-kontaktkartor i cancerceller jämfört med matchade normala celler, in för var och en av CTCF-bindningsställena (Additional file 1: Fig. S14b). Regioner med minskade interaktioner (log2 FC < -1, genomsnittlig log2-interaktion > 0) med cancerspecifika förlorade CTCF-bindningsställen och regioner med ökade interaktioner (log2 FC > 1, genomsnittlig log2-interaktion > 0) med cancerspecifika vunna CTCF-bindningsställen användes för analys av anrikning av transkriptionsfaktorer (TF) nedströms.
Transkriptionsfaktorberikningsanalys
En reviderad version av BART-algoritmen användes för TF-berikningsanalys. I korthet har en samling av union DNase I hypersensitiva platser (UDHS) tidigare kurerats som en repertoar av alla kandidat cis-regulatoriska element i det mänskliga genomet, och 7032 ChIP-seq-dataset samlades in för 883 TF:er , där varje TF har ett eller flera ChIP-seq-dataset från flera celltyper eller förhållanden. En binär profil genererades för varje TF på UDHS som anger om TF har minst en topp från någon av sina ChIP-seq-datasetterna som är lokaliserad inom varje UDHS. Analys av bindningsberikning tillämpades för varje TF genom att jämföra TF-bindningen på en delmängd av UDHS som överlappar de utvalda genomiska regionerna med TF-bindningen på UDHS. p-värdet erhölls med hjälp av Fisher’s exakta test med två svansar.