Procedimiento experimental
Xenoinjerto en pacientes y cultivo celular
Las líneas celulares humanas de T-ALL incluyen CUTLL1 (regalo de Adolfo Ferrando, Universidad de Columbia) y JURKAT (American Type Culture Collection (ATCC), Manassas, VA, #CCL-119) . Las células se cultivaron en medio RPMI1640 con l-glutamina y 25 mM de HEPES (Corning) complementado con un 10% de suero bovino fetal inactivado por calor (Sigma-Aldrich), 10 U/mL de penicilina-estreptomicina (Gibco) y 1× glutaMAX (Gibco) en una incubadora humidificada a 37 °C y 5% de CO2. Las células se analizan periódicamente para detectar la presencia de micoplasma utilizando el kit de detección de micoplasma Lonza Walkersville MycoAlert (última prueba en enero de 2020). Las líneas celulares se mantienen en cultivo durante un máximo de 20 pases y se autentifican mediante el perfil de repeticiones en tándem cortas (JURKAT) o mediante PCR para detectar la translocación TCRb-NOTCH1 (TCRBJ2S4CUTLL1F:5′-GGACCCGCTCAGTGCT-3′, NOTCH1CUTTL1R:5′-TCCCGCTCCAAAATAAGG-3′). La última autentificación celular se realizó en febrero de 2020. Las células T CD4+ humanas se compraron a AllCells. Las muestras humanas primarias se recogieron con el consentimiento informado y se analizaron bajo la supervisión de la Junta de Revisión Institucional de la Universidad de Padua, la Associazone Italiana di Ematologia e Oncologia Pediatrica, y los ensayos clínicos pediátricos ALL 2000/2006 de Berlín-Frankfurt-Münster (AIEOP-BFM). Se obtuvo el consentimiento informado para utilizar el material sobrante con fines de investigación de todos los pacientes a la entrada del ensayo de acuerdo con la Declaración de Helsinki.
Anticuerpos y reactivos
Se realizaron Western blots utilizando los siguientes anticuerpos: Actina y CTCF de Millipore Sigma (clon C4; 07-729) y NOTCH1 escindido (Val1744) de Cell Signaling Technology (4147). ChIP-seq se realizó utilizando los siguientes anticuerpos: CTCF de Millipore Sigma (07-729), H3K27Ac (8173S), y H3K27me3 (9733S) de Cell Signaling Technology, y H3K4me1 (07-473) de Millipore.
Hi-C in situ
El Hi-C in situ se realizó en células T CD4+, Jurkat, CUTLL1, y xenoinjertos de pacientes como se ha descrito previamente. En resumen, las células fueron reticuladas con formaldehído al 1% durante 10 minutos a temperatura ambiente. Por reacción Hi-C, se lisaron 5 millones de células y se permeabilizaron los núcleos. El ADN se digirió con MboI de New England Biolabs (R0147M). Los fragmentos digeridos se marcaron con d-ATP biotinilado de Jena Bioscience (NU-835-BIO14-S) y se ligaron. Tras el tratamiento con RNasa y el tratamiento con Proteinasa K para revertir los enlaces cruzados, los núcleos se sonicaron utilizando un Covaris E220 para producir una longitud de fragmento media de 400 pb. Se utilizaron perlas de estreptavidina de Thermo Fisher Scientific (65001) para arrastrar los fragmentos marcados con biotina. Tras la purificación y el aislamiento del ADN, se prepararon bibliotecas finales utilizando el NEBNext® Ultra™ II DNA Library Prep Kit for Illumina® y se secuenciaron mediante secuenciación de extremos emparejados con una longitud de lectura de 150 pb en un Illumina HiSeq 2500 para producir una media de 400 millones de lecturas por muestra.
Perfiles ChIP-seq
Las células T CD4+, Jurkat, CUTLL1 y los xenoinjertos de los pacientes fueron reticulados con formaldehído al 1% y suero bovino fetal al 1% en PBS durante 10 minutos a temperatura ambiente. La reacción se apagó con glicina 0,2 M a temperatura ambiente durante 5 minutos. A continuación, las células se lavaron con PBS y se pelletearon.
Para los ChIPs de CTCF, la inmunoprecipitación se realizó basándose en un protocolo descrito previamente. Un pellet que contenía 50 millones de células se lisó con 5 mL de tampón de lisis (50 mM de HEPES-KOH, pH 7,5, 140 mM de NaCl, 1 mM de EDTA, 10% de glicerol, 0,5% de NP-40, 0,25% de Tritón X-100) durante 10 min a 4 °C. Los núcleos se pelletearon a 1350×g durante 7 minutos y se resuspendieron en 10 mM Tris pH 8, 1 mM EDTA y 0,1% SDS. La cromatina se cortó con un sistema Covaris E220 hasta alcanzar una longitud media de fragmento de 400 pb y se centrifugó a 15.000 rpm durante 10 minutos para eliminar la cromatina insoluble y los restos. El sobrenadante se incubó con 20 μL de Dynabeads Protein G durante 30 minutos antes de desechar las perlas. El 1% del volumen total se guardó como entrada y el resto se incubó con el anticuerpo anti-CTCF durante la noche. En total, se añadieron 100 μL de Dynabeads Protein G durante 2 h. Los fragmentos unidos se lavaron dos veces con 1 mL de tampón de baja salinidad (20 mM Tris-HCl pH 8,0, 150 mM NaCl, 2 mM EDTA, 1% p/v Triton X-100 y 0.1% p/v SDS), una con tampón de alta salinidad (20 mM Tris-HCl pH 8,0, 500 mM NaCl, 2 mM EDTA, 1% p/v Triton X-100, y 0,1% p/v SDS), una con tampón de cloruro de litio (10 mM Tris-HCl pH 8.0, 250 mM LiCl, 1 mM EDTA, 1% p/v NP-40, y 1% p/v de ácido desoxicólico), y dos veces con TE (10 mM Tris pH 8, 1 mM EDTA).
Para los ChIP de histonas, las células se lisaron en 375 μL de tampón de incubación de núcleos (15 mM Tris pH 7.5, 60 mM KCl, 150 mM NaCl, 15 mM MgCl2, 1 mM CaCl2, 250 mM de sacarosa, 0,3% de NP-40, 1 mM de NaV, 1 mM de NaF y 1 pastilla de inhibidor de proteasas sin EDTA (Roche)/10 mL en H2O) durante 10 min en hielo. Los núcleos se lavaron una vez con tampón de digestión (10 mM NaCl, 10 mM Tris pH 7,5, 3 mM MgCl2, 1 mM CaCl2, 1 mM NaV, 1 mM NaF y 1 pastilla inhibidora de proteasas sin EDTA (Roche)/10 mL en H2O) y se resuspendieron en 57-μL de tampón de digestión que contenía 4,5 unidades de MNasa (USB) durante 1 h a 37 °C. La actividad de la MNasa se apagó durante 10 minutos en hielo tras la adición de EDTA hasta una concentración final de 20 mM. Los núcleos se pelletearon y se resuspendieron en 300-μL de tampón de lisis de núcleos (50 mM Tris-HCl pH 8,0, 10 mM EDTA pH 8,0, 1% SDS, 1 mM NaV, 1 mM NaF y 1 pastilla inhibidora de proteasas sin EDTA (Roche)/10 mL en H2O) antes de la sonicación con un Bioruptor Pico (Diagenode) durante 5 min (30 s encendido, 30 s apagado). El lisado se centrifugó a velocidad máxima durante 5 minutos para eliminar los restos. Se añadieron al sobrenadante nueve volúmenes de tampón de dilución IP (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 y 1 pastilla de inhibidor de proteasa sin EDTA (Roche)/10 mL en H2O). En total, se añadieron 50 μL de Dynabeads Protein G y la muestra se incubó a 4 °C durante 30 minutos, rotando. El 1% de la muestra se mantuvo como entrada, y la muestra restante se dividió en 3 tubos. En total, se añadieron 50 μL de Dynabeads Protein G conjugadas con 15 μL del anticuerpo apropiado a cada tubo antes de la incubación de una noche a 4 °C, rotando. Los complejos unidos a las microesferas se lavaron durante 5 minutos cada uno en 1 mL de tampón de baja salinidad, tampón de alta salinidad, tampón LiCl y dos veces con TE.
Para eluir los complejos unidos a las microesferas, éstas se resuspendieron en 50 μL de tampón de elución (100 mM de NaHCO3, 1% p/v de SDS) y se incubaron a 65 °C durante 15 min, agitando a 1000 RPM en un termomixer (Thermo Scientific). La elución se repitió una segunda vez y, a continuación, se añadieron 100 μL de tampón RNasa (12 μL de NaCl 5 M, 0,2 μL de RNasa 30 mg/mL y 88 μL de TE) a cada muestra de ChIP y de entrada. Las muestras se incubaron a 37 °C durante 20 minutos, seguidas de la adición de 100 μL de tampón de proteinasa K (2,5 μL de 20 mg/mL de proteinasa K, 5 μL de SDS al 20% y 92,5 μL de TE) durante la noche a 65 °C. Se añadió un volumen igual de solución de fenol:cloroformo y se mezcló bien. La mezcla se transfirió a tubos MaXtract de alta densidad (Qiagen) y se centrifugó durante 8 minutos a 15.000 rpm. La fase superior se transfirió a nuevos tubos y se mezcló con 1,5 μL de 20 mg/mL de glucógeno, 30 μL de acetato de sodio 3M y 800 μL de etanol. Las muestras se incubaron a – 80 °C hasta que se congelaron y luego se centrifugaron a 15.000 rpm durante 30 min a 4 °C. Se eliminó el sobrenadante y los pellets se lavaron en 800 μL de etanol al 70% helado y se centrifugaron durante 10 min a 4 °C a 15.000 rpm. Tras una cuidadosa eliminación del etanol, los pellets se secaron al aire y se resuspendieron en 30 μL de Tris 10 mM a pH 8.
Las IP y el ADN de entrada se cuantificaron entonces utilizando un fluorímetro Qubit 3.0. Las bibliotecas se prepararon utilizando el KAPA HyperPrep Kit (KK8505) y se secuenciaron con un Illumina NextSeq 500 hasta una profundidad media de 28 millones de lecturas por muestra.
Perfiles de ARN-seq
El ARN se aisló de 3 millones de células por muestra utilizando el Bio-Rad Aurum™ Total RNA Mini Kit y se cuantificó con el Agilent RNA 6000 Nano Kit con el Agilent Bioanalyzer. Las bibliotecas se prepararon mediante el agotamiento del ARNr utilizando el Illumina TruSeq® Stranded mRNA Library Prep Kit para una baja concentración de muestra inicial y se secuenciaron mediante secuenciación de extremo único en un Illumina NextSeq 500 hasta una profundidad media de 18 millones de lecturas por muestra.
Perfil de metilación del ADN
El ADN genómico se aisló utilizando el AllPrep DNA/RNA Micro Kit (Qiagen). Para evaluar el estado de metilación del ADN en todo el genoma, realizamos mRRBS . Tras la cuantificación fluorométrica utilizando un instrumento Qubit 3.0, digerimos el ADN genómico con la enzima de restricción MspI (New England Biolabs) y seleccionamos por tamaño fragmentos de aproximadamente 100-250 pares de bases de longitud utilizando perlas de inmovilización reversible en fase sólida (SPRI) (MagBio Genomics). El ADN resultante se sometió a la conversión de bisulfito utilizando el kit EZ DNA Methylation-Lightning (Zymo Research). Creamos bibliotecas a partir de ADN monocatenario convertido por bisulfito utilizando el Pico Methyl-Seq Library Prep Kit (Zymo Research), que luego se agruparon para la secuenciación en un instrumento Illumina NextSeq 500 utilizando el kit de reactivos NextSeq 500/550 V2 High Output (1 × 75 ciclos) hasta una profundidad de lectura mínima de 50 millones de lecturas por muestra.
Secuenciación del genoma completo
Tres millones de células procedentes de líneas celulares o de muestras de pacientes fueron peleadas y resuspendidas en 1 mL de Cell Lysis Solution (Qiagen) mezclada con 500 μg de RNasa A. La reacción de lisis se llevó a cabo a 37 °C durante 15 min. En total, se añadieron 333 μL de solución de precipitación de proteínas (Qiagen) a cada muestra, que se agitó en vórtex y luego se centrifugó a 2000×g durante 10 min. El sobrenadante se mezcló con 1 mL de isopropanol hasta que las hebras de ADN se precipitaron de la solución. Tras desechar el sobrenadante, el pellet de ADN se lavó con 1 mL de etanol al 70% y se centrifugó a 2000×g durante 1 min. A continuación se vertió el etanol y el pellet se secó al aire durante 15 minutos antes de resuspenderlo en 50 a 100 μL de solución de hidratación de ADN (Qiagen). El ADN se secuenció con la secuenciación Illumina de extremos pareados con una cobertura de 30×.
Inmunoprecipitación
Se pelleteó un total de 100 millones de células para cada reacción de inmunoprecipitación y se incubó en el tampón A (10 mM HEPES pH 8,0, 1,5 mM MgCl2, 10 mM KCl, 0,5 mM DTT) durante 10 min en hielo. A continuación, las células se lisaron tras 12 golpes con una trituradora de tejidos de 7 mL con mortero suelto (Wheaton, 357542) y se centrifugaron a 2000 rpm durante 7 min. Los pellets nucleares se resuspendieron en 5 volúmenes de tampón TENT (50 mM Tris pH 7,5, 5 mM EDTA, 150 mM NaCl, 1% Triton X-100, 5 mM MgCl2) y se trataron con benzonasa durante 30 min antes de realizar 5 pases por una jeringa de 25 g × 5/8 pulg. La fracción insoluble se eliminó tras centrifugación a 2000 rpm durante 7 min y se incubó durante la noche con Dynabeads Protein G hibridada con anticuerpos. Se extrajo un total de 2 millones de células para la entrada. Los lisados de perlas y núcleos se lavaron 6 veces con tampón TENT y luego se eluyeron en glicina 0,1 M de pH 2,5 con 100 mM de Tris de pH 8,0 previamente. Se añadió tampón de muestra NuPAGE LDS a los eluidos y a las entradas, que se incubaron a continuación a 70 °C durante 15 min antes del análisis por western blot.
Recogida de datos públicos
Los datos públicos de CTCF ChIP-seq se recogieron de Cistrome Data Browser (para archivos de picos) y NCBI GEO (para archivos fastq, Archivo adicional 2: Tabla S1). Los datos de ChIP-seq de modificación de histonas se recogieron de NCBI GEO y ENCODE (para archivos bam). Los datos públicos de RNA-seq en múltiples tipos de células se recogieron de ENCODE (para archivos fastq). Los datos del perfil de metilación del ADN se recogieron de ENCODE (para archivos bedMethyl) y NCBI GEO. Los datos de Hi-C se recogieron de NCBI GEO y ENCODE (para archivos fastq). Los datos de ATAC-seq se recogieron de NCBI GEO (para los archivos fastq). Los datos de secuenciación del genoma completo de las muestras BRCA, COAD, LUAD y PRAD se recogieron del Portal de Datos del Consorcio Internacional del Genoma del Cáncer (ICGC). La información detallada, incluyendo los ID de acceso de todos los conjuntos de datos públicos recogidos en este trabajo, se puede encontrar en el archivo adicional 6: Tabla S5.
Procesamiento de datos
Análisis de datos ChIP-seq
La alineación de secuencias para los datos ChIP-seq en archivos fastq se llevó a cabo utilizando la misma tubería de análisis estándar que se utiliza en Cistrome DB , para la consistencia y reproducibilidad. Todos los datos de secuencia de alineación genómica se realizaron utilizando la tubería Chilin con los parámetros por defecto ($ chilin simple -p estrecha -s hg38 –threads 8 -t IN.fq -i PRENAME -o OUTDIR). Brevemente, las lecturas de la secuencia se alinearon con el genoma de referencia humano (GRCH38/hg38) utilizando 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). A continuación, los archivos Sam se convirtieron en archivos bam utilizando samtools ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). Para los conjuntos de datos CTCF ChIP-seq, se utilizó MACS2 para llamar a los picos bajo el umbral FDR de 0,01 ($ macs2 callpeak –SPMR -B -q 0,01 –keep-dup 1 -g hs -t PRENAME.bam -n PRENAME –outidr OUTDIR). Se retuvieron los picos con un enriquecimiento de pliegues de al menos 4. Los archivos Bigwiggle se generaron utilizando las herramientas BEDTools y UCSC ($ 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). Por último, sólo se incluyeron en el análisis integrador posterior las muestras de CTCF ChIP-seq que tenían al menos 2000 picos.
Análisis de datos de ATAC-seq
Se utilizó Trim Galore para recortar las lecturas de secuenciación en bruto ($ trim_galore –nextera –phred33 –fastqc –paired R1.fq R2.fq -o OUTDIR). Las lecturas se alinearon con el genoma de referencia humano (GRCH38/hg38) utilizando Bowtie2 ($ bowtie2 -p 10 -X 2000 -x INDEX -1 R1.fq -2 R2.fq -S PRENAME.sam). A continuación, los archivos Sam se convirtieron en archivos bam utilizando samtools ($ samtools view -bS -q 1 -@ 8 PRENAME.sam > PRENAME.bam). Se utilizó Bedtools para convertir los archivos bam en formato bed ($ bamToBed -i PRENAME.bam -bedpe > PRENAME_PE.bed). Las lecturas asignadas al ADN mitocondrial se descartaron del análisis posterior.
Análisis de datos de ARN-seq
Los conjuntos de datos de ARN-seq se procesaron con Salmon ($ salmon quant –gcBias -i INDEX -l A -p 8 {-1 R1.fq -2 R2.fq| -r IN.fq} -o OUTDIR). El índice del transcriptoma se construyó sobre el genoma humano de referencia (GRCH38/hg38). Las estimaciones de abundancia a nivel de transcripción se resumieron a nivel de genes utilizando el paquete «tximport» para el análisis de expresión diferencial. DESeq2 se utilizó para identificar los genes expresados diferencialmente, y los diferentes umbrales utilizados en diferentes análisis se enumeran correspondientemente en el manuscrito.
Análisis de datos Hi-C
Los datos Hi-C se procesaron utilizando HiC-Pro ($ HiC-Pro -i INDIR -o OUTDIR -c CONFIG -p). Se generaron mapas de contacto con una resolución de 5 kb. Los datos de la matriz en bruto se normalizaron utilizando el enfoque descrito en Normalización de las interacciones de la cromatina.
Análisis de los datos de metilación del ADN
Los datos de metilación del ADN (para las líneas celulares de T-ALL y los pacientes de T-ALL) se desmultiplexaron con bcl2fastq seguido de un recorte de 10 pares de bases del extremo 5′ para eliminar las secuencias de cebadores y adaptadores utilizando TrimGalore . La alineación de la secuencia con el genoma de referencia GRCh38/hg38 y las llamadas de metilación se realizaron con Bismark ($ bismark –multicore 8 –bowtie2 -q -N 1 INDEX INFILE.fq). Los archivos de cobertura (recuentos) para las citosinas en el contexto CpG se generaron con Bismark ($ bismark_methylation_extractor –multicore 8 –comprehensive –bedGraph INFILE_bismark_bt2.bam).
Análisis de datos de secuenciación del genoma completo
Se identificaron mutaciones para dos líneas celulares de T-ALL (Jurkat y CUTLL1) y dos muestras de pacientes de T-ALL a partir de los datos de secuenciación del genoma completo. Se alinearon las secuencias de hilo corto de Illumina con el genoma humano de referencia (GRCH38/hg38) utilizando BWA mem. Utilizamos SAMBlaster para identificar los pares discordantes, dividir las lecturas y marcar los posibles duplicados de PCR. Utilizamos SAMBAMBA para convertir los SAM alineados al formato BAM, y se utilizó samtools para ordenar los alineados y crear un archivo BAM correspondiente a cada muestra.
Utilizamos VarDict para identificar las variantes que se solapaban con los sitios de unión CTCF. Utilizamos todos los parámetros por defecto, excepto «-f 0,1», que se utilizó para identificar las variantes que fueron apoyadas por más del 10% de las lecturas en ese lugar. Anotamos las variantes utilizando Variant Effect Predictor (VEP) y utilizamos scripts personalizados para identificar las variantes que influyen en la unión al TF.
Utilizamos de nuevo VarDict para identificar las variantes en los genes CTCF y NOTCH1 para las cuatro muestras. Utilizamos todos los parámetros por defecto excepto «-f 0.1» que se utilizó para identificar las variantes que fueron apoyadas por más del 10% de las lecturas en esa ubicación. Anotamos las variantes utilizando Variant Effect Predictor (VEP) , y luego lo filtramos para identificar las mutaciones que o bien (a) no se veían en más del 1% de cualquier población humana normal, o (b) tenían una puntuación CADD de deletéreo > 20, o (c) estaba presente en la base de datos COSMIC.
Modelado integral y análisis estadístico
Identificación del repertorio de unión de CTCF en el genoma humano
Para la ChIP-seq de CTCF, recogimos un total de 793 conjuntos de datos, incluyendo 787 conjuntos de datos públicos y 6 conjuntos de datos generados por nosotros (Archivo adicional 2: Tabla S1). En total, se utilizaron en este estudio 771 conjuntos de datos CTCF ChIP-seq con picos superiores a 2000. Cada conjunto de datos puede producir picos CTCF identificados por MACS2 en un rango entre 2050 y 198.021, con una mediana de 46.451 y un total de 36.873.077 picos (archivo adicional 1: Fig. S1a). La distribución de las longitudes de los intervalos entre cimas de picos CTCF adyacentes de todos los 36.873.077 picos de los 771 conjuntos de datos tiene un punto de inflexión a ~ 150 pb (archivo adicional 1: Fig. S1c) que indica el límite entre el mismo sitio de unión y diferentes sitios de unión . Por lo tanto, utilizamos 150 pb como límite para fusionar los picos de CTCF. En la práctica, extendimos ± 75 pb desde la cima de cada pico para generar una región de 150 pb centrada en la cima para representar cada pico y fusionamos todas las regiones de picos superpuestos para generar un conjunto de unión de sitios de unión de CTCF, que contiene 688.429 sitios no superpuestos. A cada sitio de unión se le asignó una puntuación de ocupación de CTCF, definida como el recuento de conjuntos de datos ChIP-seq que muestran un pico dentro del sitio. En consecuencia, definimos la frecuencia de ocupación como la relación de la puntuación de ocupación sobre el número total de conjuntos de datos ChIP-seq de CTCF. Para garantizar aún más la solidez de los sitios de unión CTCF identificados, seleccionamos 285.467 sitios de alta confianza con una puntuación de ocupación ≥ 3 para los análisis posteriores. Los motivos CTCF dentro de los sitios de unión se buscaron mediante FIMO con la matriz Jaspar (ID: MA0139.1), con un umbral de valor p de 1e-4. Se retuvo un motivo con el menor valor p para cada sitio de unión a CTCF.
Identificación de sitios de unión a CTCF constitutivos
La distribución de las puntuaciones de ocupación de todos los 285.467 sitios de unión a CTCF (Archivo adicional 1: Fig. S1d, curva azul) muestra que la mayoría de los sitios de unión a CTCF ocurren sólo en unos pocos conjuntos de datos, y el número de sitios de unión disminuye con el aumento de la puntuación de ocupación cuando ésta es pequeña. Sin embargo, hay sitios de unión de CTCF que están muy conservados en casi todos los conjuntos de datos (por ejemplo, sitios de unión con una puntuación de ocupación superior a 600). Utilizamos una función de ley de potencia para ajustar la curva de distribución (azul) que se muestra en el archivo adicional 1: Fig. S1d para determinar el límite de los sitios CTCF constitutivos. Denotamos Oi como el número de sitios de unión CTCF observados con puntuación de ocupación igual a i, y Ei como el número de sitios CTCF esperados con puntuación de ocupación igual a i. El ajuste de la ley de potencia a los datos Oi puede describirse como (archivo adicional 1: Fig. S1d, verde):
Definimos el punto de corte A para los sitios de unión constitutivos de CTCF como:
En otras palabras, el total de sitios CTCF observados con una puntuación de ocupación mayor que A debería ser 6 veces más de lo esperado. Entonces determinamos A = 615, y utilizamos un corte de frecuencia de ocupación del 80% para definir 22.097 sitios de unión a CTCF constitutivos, que corresponde a la puntuación de ocupación ≥ 616 en todos los 771 conjuntos de datos ChIP-seq de CTCF.
Identificación de sitios de unión a CTCF ganados/perdidos específicos del cáncer
Utilizamos los siguientes 2 criterios para identificar sitios de unión a CTCF perdidos específicos del cáncer: (1) El sitio de unión a CTCF debe tener una frecuencia de ocupación más baja para los conjuntos de datos de ese tipo de cáncer en comparación con la frecuencia de ocupación para todos los conjuntos de datos y (2) el nivel de unión a CTCF (cuantificado como recuentos de lecturas ChIP-seq normalizados) en el sitio es menor en los conjuntos de datos de cáncer que en otros conjuntos de datos. Para los sitios CTCF ganados, utilizamos el conjunto de criterios viceversa. Brevemente, para cada sitio de unión a CTCF en cada tipo de cáncer, se calculó la puntuación de ocupación en los conjuntos de datos de cáncer junto con su puntuación de ocupación en los 771 conjuntos de datos. Los niveles de unión a CTCF se obtuvieron a partir de una matriz de recuento de lecturas normalizada en la que primero se calcularon los recuentos de lecturas de ChIP-seq (RPKM) para los sitios de unión a CTCF en todos los conjuntos de datos y luego se realizó una normalización por cuantiles. Se utilizó la prueba t de Student no apareada de dos colas para cuantificar la diferencia de los niveles de unión entre los diferentes grupos de conjuntos de datos, y el valor p se ajustó entonces utilizando el procedimiento de Benjamini-Hochberg . Además, se compararon las puntuaciones de ocupación de la unión y los niveles de unión entre los conjuntos de datos de cáncer y los conjuntos de datos de los tejidos o tipos de células normales emparejados, con el fin de tener en cuenta el posible factor de confusión de la especificidad del tejido en lugar de la especificidad del cáncer. A continuación se describen los criterios detallados para identificar los sitios de unión de CTCF específicos del cáncer:
-
Sitios de unión de CTCF perdidos específicos del cáncer: (1) frecuencia de ocupación ≤ 0,2 en conjuntos de datos de cáncer; (2) frecuencia de ocupación ≥ 0,7 en 771 conjuntos de datos; (3) frecuencia de ocupación ≥ 0.5 (con puntuación de ocupación ≥ 2) en conjuntos de datos de tejido normal emparejados; (4) los niveles de CTCF son más bajos en el cáncer en comparación con todos los demás conjuntos de datos (puntuación estadística < 0), (5) los niveles de CTCF son más bajos en el cáncer en comparación con los conjuntos de datos de tejido normal emparejados (puntuación estadística < 0), (6) las señales de unión de CTCF promediadas (RPKM) < 5 en los conjuntos de datos de cáncer.
-
Sitios de unión a CTCF ganados específicamente para el cáncer: (1) frecuencia de ocupación ≥ 0,5 (con puntuación de ocupación ≥ 2) en conjuntos de datos de cáncer, (2) frecuencia de ocupación ≤ 0,2 en 771 conjuntos de datos, (3) puntuación de ocupación = 0 en conjuntos de datos de tejido normal emparejados, (4) los niveles de CTCF son significativamente mayores en cáncer en comparación con todos los demás conjuntos de datos (FDR ≤ 0.01), (5) los niveles de unión de CTCF son significativamente más altos en el cáncer en comparación con los conjuntos de datos de tejido normal emparejados (FDR ≤ 0,01), (6) las señales de unión de CTCF promediadas (RPKM) > 2 en los conjuntos de datos de cáncer.
Los sitios de unión a CTCF específicos ganados y perdidos para cada tipo de cáncer se muestran en el archivo adicional 4: Tabla S3.
Cuantificación de la accesibilidad diferencial a la cromatina
Utilizamos los datos procesados de la Ref. que incluyen una matriz de recuentos de inserción ATAC-seq normalizados dentro del conjunto de picos de cáncer TCGA para evaluar la accesibilidad diferencial a la cromatina alrededor de los sitios de unión a CTCF. Para cada tipo de cáncer entre BRCA, CRC, LUAD y PRAD, los picos de ATAC-seq pan-cáncer que se solapan con los sitios de unión a CTCF perdidos o ganados específicos del tipo de cáncer se utilizaron para los análisis posteriores. La puntuación diferencial de ATAC-seq para cada pico se cuantificó como el cambio de pliegues de la media de los recuentos de inserción de ATAC-seq normalizados de las muestras de pacientes en el tipo de cáncer correspondiente frente a las de pacientes en otros tipos de cáncer, y la puntuación diferencial de ATAC-seq se asignó entonces al pico que se solapaba con el sitio de unión a CTCF.
Para mantener la coherencia, aplicamos el mismo enfoque utilizado para los datos de ATAC-seq del TCGA para analizar los datos de ATAC-seq recogidos de la línea celular T-ALL Jurkat y las células T CD4+ normales. Se generó una matriz de datos utilizando los recuentos de lecturas brutas de ATAC-seq en los sitios de unión de CTCF para todos los conjuntos de datos de Jurkat y células T. Se aplicó la normalización cuantílica en la matriz a escala log2 (pseudoconteo = 5). La puntuación diferencial ATAC-seq se midió como el cambio de pliegue de los recuentos ATAC-seq normalizados promediados entre los conjuntos de datos de Jurkat frente a los de células T CD4+ en cada sitio de unión CTCF.
Normalización de las interacciones de la cromatina
Dado un mapa de contactos Hi-C A = {aij}, la puntuación aij refleja las lecturas mapeadas entre dos regiones genómicas i y j. Supongamos que el tamaño del bin es de 5 kb, las regiones i y j tendrán una distancia genómica de ∣i – j ∣ × 5kb. Dado que la probabilidad de contacto entre dos bin disminuye con el aumento de la distancia genómica , normalizamos el mapa de contactos de la siguiente manera: para cualquier distancia genómica dada dk = k × 5kb, cuantificamos un factor de normalización \( {\overline{S}_{d_k} \) como las interacciones promediadas entre todos los pares de bin con la misma distancia genómica dk en un mismo cromosoma, por ejemplo \( {\overline{S}_{d_k}=\left({\sum}_{j-i=k}{a}_{ij}\right)/n \n), donde n es el número total de pares de contenedores con distancia dk. La puntuación de la interacción aij entre dos recipientes con distancia dk se normalizó mediante \( {\overline{S}_{d_k} \) como \( {a}_{ij}^{prime }={a}_{ij}/{overline{S}_{d_k} \). Utilizando este enfoque, normalizamos la matriz A en \( A^{{prime }={left{{a}_{ij}^{{prime}{right}}) dentro de cada cromosoma.
Detección de interacciones diferenciales de la cromatina
Denotamos los mapas de contacto Hi-C normalizados en el conjunto de datos de cáncer y el conjunto de datos normal como C = {cij} y N = {nij}, respectivamente. Para un sitio de unión CTCF dado x (con la coordenada xc) y una distancia genómica predefinida L, las interacciones de la cromatina entre x y sus bines cercanos no solapados de 5 kb con distancia genómica hasta L se recogen de C y N respectivamente. Específicamente, las puntuaciones de las interacciones entre x y sus fragmentos cercanos de 5 kb en C se recogen como IC = {cij} mientras que i o j es igual a ⌊xc/5kb⌋, y 0 < (j – i) × 5kb ≤ L. Del mismo modo, las puntuaciones de interacción entre x y su cercano 5-kb bins en N se recogieron como IN = {nij}. A continuación, se aplicó una prueba t de Student emparejada de dos colas sobre IC e IN para cuantificar la interacción diferencial entre el cáncer y las células normales que rodean el sitio de unión a CTCF x.
Asociación de la unión de CTCF con la expresión génica
En total, se seleccionaron 54 tipos celulares para los que tanto los datos ChIP-seq de CTCF como los datos RNA-seq están disponibles públicamente (archivo adicional 6: Tabla S5) para investigar la asociación entre la unión de CTCF y la expresión génica para cada par CTCF-gen en el mismo cromosoma. Para obtener el nivel de unión de CTCF, se generó una matriz de recuento de lecturas utilizando lecturas por kilobase por millón (RPKM) en los sitios de unión de CTCF de los datos de ChIP-seq. La matriz de recuento de lecturas se escaló con la raíz cuadrada de RPKM seguida de una normalización cuantitativa. El nivel de expresión génica se midió para cada gen utilizando la raíz cuadrada de los transcritos por millón (TPM) de los datos de RNA-seq. Para cada par CTCF-gen, cuantificamos la asociación entre el sitio CTCF y el gen en los 54 tipos de células utilizando el coeficiente de correlación R entre el nivel de unión CTCF normalizado y la expresión génica (Fig. 3a). Los pares CTCF-gen se consideraron «altamente correlacionados» con un R2 superior a 0,25, es decir coeficiente de correlación superior a 0,5 o inferior a – 0,5, y los pares CTCF-genes altamente correlacionados contribuyen al 1,3% de todos los pares CTCF-genes (archivo adicional 1: Fig. S8a).
Identificación de los dominios de cromatina limitados por CTCF constitutivo
Para cada sitio de unión a CTCF, definimos su dominio de cromatina asociado como la región genómica que (1) incluye este sitio de unión a CTCF específico, (2) está limitada por un par de sitios de unión a CTCF constitutivo con motivos de orientaciones opuestas, y (3) ocupa un mínimo de 100 kb y un máximo de 1 MB de región a cada lado del sitio de unión a CTCF. La Figura 3b contiene un esquema de cómo se definieron los dominios de cromatina constitutivos de CTCF.
Detección de los cambios de metilación del ADN que rodean los sitios de unión a CTCF
Los cambios de metilación del ADN se detectaron dentro de una región de 300 pb centrada en cada sitio de unión a CTCF. Se retuvieron las regiones con al menos 3 CpGs cubiertas por al menos 5 lecturas (≥ 5×) tanto en las líneas celulares de cáncer como en los tejidos normales correspondientes. Una región de 300 pb se detectó como diferencialmente metilada si el promedio de los niveles de metilación diferencial de todos los CpGs (≥ 5×) dentro de esta región era superior al 20%.
Detección de la tasa de mutación y la puntuación del motivo diferencial
Para cada sitio de unión a CTCF, el recuento de mutaciones en bruto se calculó como la aparición de eventos de mutación en todas las muestras/pacientes en cada par de bases individuales dentro de una región de 400 pb centrada en el sitio de unión a CTCF. La tasa de mutación para un grupo de sitios de unión a CTCF se calculó como el recuento promedio de mutaciones sobre el número de sitios de unión a CTCF para cada par de bases dentro de la región de 400 pb.
La puntuación del motivo se midió mediante la puntuación de la matriz de peso de la posición de CTCF (Jaspar , Matrix ID: MA0139.1) a una secuencia de ADN de 19 pb centrada en el motivo de CTCF o en el sitio de unión a CTCF utilizando cocientes de probabilidad logarítmica (con la frecuencia de nucleótidos de fondo como para A,C,G,T). La puntuación diferencial de los motivos se calculó comparando las puntuaciones de los motivos para las secuencias de referencia y las mutadas.
Análisis de motivos de la secuencia de ADN
El análisis de enriquecimiento de motivos de la secuencia de ADN se realizó utilizando MDSeqPos (versión 1.0.0) en Cistrome con los parámetros por defecto (-cisrome -Homo Sapien o Mus musculus). Los análisis de motivos de novo se realizaron utilizando HOMER (versión 4.10) con el módulo findmotifs.pl y MEME (versión 5.1.).1) con los siguientes parámetros: meme -dna -mod zoops -maxw 20 -evt -0,01.
Identificación de las regiones de interacción diferencial intradominio de CTCF
Para un conjunto dado de sitios de unión de CTCF, se recogieron los cambios de interacción de la cromatina entre un sitio de CTCF y cada uno de sus bines no solapados intradominio, medidos a partir de mapas de contacto Hi-C normalizados en células cancerosas sobre células normales emparejadas, para cada uno de los sitios de unión de CTCF (archivo adicional 1: Fig. S14b). Las regiones con interacciones disminuidas (log2 FC < -1, interacción log2 promediada > 0) con sitios de unión CTCF perdidos específicos del cáncer, y las regiones con interacciones aumentadas (log2 FC > 1, interacción log2 promediada > 0) con sitios de unión CTCF ganados específicos del cáncer se utilizaron para el análisis de enriquecimiento del factor de transcripción (TF) aguas abajo.
Análisis de enriquecimiento de factores de transcripción
Se utilizó una versión revisada del algoritmo BART para el análisis de enriquecimiento de TF. Brevemente, se curó previamente una colección de sitios hipersensibles a la DNasa I de la unión (UDHS) como un repertorio de todos los elementos cis-reguladores candidatos en el genoma humano, y se recogieron 7032 conjuntos de datos ChIP-seq para 883 TFs , teniendo cada TF uno o más conjuntos de datos ChIP-seq de múltiples tipos de células o condiciones. Se generó un perfil binario para cada TF en UDHS indicando si el TF tiene al menos un pico de cualquiera de sus conjuntos de datos ChIP-seq localizados dentro de cada uno de los UDHS. Se aplicó un análisis de enriquecimiento de la unión para cada TF comparando la unión del TF en un subconjunto de UDHS que se solapaba con las regiones genómicas seleccionadas frente a la unión del TF en los UDHS. El valor p se obtuvo utilizando la prueba exacta de Fisher de dos colas.