¿Por qué los llamadores de SV pierden las CNV grandes?

En la entrega anterior, analizamos las SV con Manta, Delly y GRIDSS. Sin embargo, estas herramientas pasan por alto ciertos tipos de SV en la práctica clínica. Se trata de deleciones (DEL) y duplicaciones (DUP) muy grandes, especialmente aquellas de más de 100 kbp.

La razón es que estas herramientas utilizan como firma principal las lecturas discordantes de extremo emparejado (paired-end) y las lecturas divididas (split reads). Cuando hay una gran DEL, el número de lecturas involucradas en cada punto de ruptura es bajo, lo que resulta en una potencia estadística insuficiente. Por el contrario, la señal más distintiva que esta gran DEL genera en el genoma es una disminución de la cobertura. Por ello, las herramientas específicas para CNV utilizan la profundidad de lectura (read depth) como su firma principal.

Tipos de CNV y su relevancia práctica

• CNV germinales: Polimorfismos naturales (por ejemplo, CNV del gen del receptor de LDL). Investigación en genética de poblaciones y trastornos del desarrollo.
• CNV somáticas: Amplificaciones génicas (por ejemplo, HER2 · MYC) o deleciones (por ejemplo, PTEN · CDKN2A) comunes en el cáncer. Indicadores para la selección de terapias dirigidas.
• CNV pequeñas a nivel de exón: Deleción de algunos exones dentro de un gen en WES. Este tipo representa el 20～30% de los diagnósticos de enfermedades genéticas.

Cada tipo requiere herramientas diferentes.

Visión general del enfoque basado en la profundidad de lectura

El principio se divide en tres pasos.

1. División del genoma en bins

Se divide todo el genoma en bins de un tamaño fijo (por ejemplo, 1 kb) y se cuenta el número de lecturas en cada bin. Este conteo es la profundidad bruta (raw depth) del bin.

2. Corrección de sesgos

La profundidad bruta presenta varios sesgos:

• Sesgo GC: Los bins con extremos de contenido de GC se amplifican menos durante la preparación de la biblioteca → menor profundidad.
• Sesgo de mapeabilidad: Los bins cerca de secuencias repetitivas no permiten que las lecturas se unan de forma única → menor profundidad.
• Momento de replicación (replication timing): Sesgo del ciclo celular en el momento de la preparación de la biblioteca.

Estos sesgos se corrigen mediante regresión local o GAM, entre otros métodos.

3. Segmentación

A partir del perfil de profundidad corregido, se identifican los puntos de ruptura donde los valores de los bins adyacentes cambian abruptamente y se agrupan los intervalos con valores similares en segmentos. Estos definen las coordenadas de inicio y fin de la CNV.

Si la profundidad media por segmento es la mitad de la referencia, es una DEL heterocigota; si es 0, una DEL homocigota; si es 1.5×, una DUP heterocigota; y si es 2×, una DUP homocigota.

CNVkit: el estándar que soporta tanto WES como WGS

Herramienta de Python creada por Eric Talevich. Está especializada en CNV somáticas (clínica oncológica), pero también puede utilizarse para germinales.

Flujo de ejecución

# 1) Crear una referencia de CNV con varias muestras normales
cnvkit.py batch \
    -n normal1.bam normal2.bam normal3.bam ... \
    -f hg38.fa \
    -t exome_targets.bed \
    -a exome_antitargets.bed \
    --output-reference myref.cnn \
    -d ref_output/

# 2) Llamar CNV tumorales
cnvkit.py batch tumor1.bam tumor2.bam \
    -r myref.cnn \
    -d results/

• Lecho objetivo (-t): región diana (lecho de exones en WES).
• Lecho antidiana (-a): región fuera de la diana (utiliza lecturas off-target).
• Referencia (myref.cnn): línea base construida con muestras normales.

La idea de CNVkit: incluso las lecturas off-target se utilizan como señal de CNV. El WES tiene lecturas en regiones fuera de la diana; no usarlas sería un desperdicio.

Resultados

Relación log2 del número de copias por segmento en el archivo .cns.

• log2=0: igual que la referencia.
• log2=−1: la mitad (candidato a delección heterocigota).
• log2=+0.58: 1,5 veces (candidato a duplicación heterocigota).
• log2=+1: el doble (candidato a duplicación homocigota).

La visualización genera perfiles de número de copias por gen en PDF mediante cnvkit.py scatter y cnvkit.py diagram.

GATK CNV — Sistema de mejores prácticas

GATK4 del Broad Institute también incluye herramientas para CNV. Enfoque basado en un Panel de Normales (PoN).

Germinales

gatk GermlineCNVCaller \
    --input case.hdf5 \
    --contig-ploidy-calls ploidy_calls/ \
    --model cohort-model \
    --output cnv_output

Se necesita un modelo creado con múltiples muestras normales.

Somático

# 1) Obtener los recuentos de lecturas del par tumor-normal
gatk CollectReadCounts -I tumor.bam -L intervals.bed --format HDF5 -O tumor.hdf5
gatk CollectReadCounts -I normal.bam -L intervals.bed --format HDF5 -O normal.hdf5

# 2) Normalizar con el PoN de muestras normales
gatk DenoiseReadCounts \
    -I tumor.hdf5 \
    --count-panel-of-normals pon.hdf5 \
    --standardized-copy-ratios tumor.std.tsv \
    --denoised-copy-ratios tumor.denoised.tsv

# 3) Segmentar y llamar
gatk ModelSegments --denoised-copy-ratios tumor.denoised.tsv -O tumor.seg
gatk CallCopyRatioSegments -I tumor.seg -O tumor.called.seg

GATK CNV está integrado en el flujo de trabajo en la nube (Terra), por lo que se integra sin problemas si es un canalón estándar del equipo.

Control-FREEC · cn.mops — Alternativas

• Control-FREEC: Herramienta para detectar CNV en pares tumor-normal desarrollada por un instituto de investigación francés. También integra información de frecuencia de alelo B (BAF).
• cn.mops: Paquete de R. Realiza la llamada conjunta integrando la profundidad de lectura de múltiples muestras mediante un modelo de mezcla de Poisson. Es estable para cohortes pequeñas.

Matriz de selección de herramientas

Corrección del sesgo GC — Perspectiva práctica

Intentemos captar la esencia del sesgo GC a través de la intuición. El GC medio del genoma humano es del 40%. El porcentaje de GC en cada bin de 1 kb se distribuye entre el 20 % y el 70 %. La curva de sesgo varía según la librería.

Librerías sin PCR: Casi no hay sesgo GC. Corrección mínima. Librerías con PCR: Máximo alrededor del 40-50 % de GC, disminuyendo fuera de ese rango (curva en forma de U). Nextera XT: Disminución extrema por debajo del 20 % de GC y por encima del 70 %.

La opción --drop-low-coverage de CNVkit elimina automáticamente estos bins extremos. En librerías con un sesgo GC pronunciado, el perfil log2 queda cubierto de ruido; por tanto, se debe realizar la segmentación después de la corrección.

Segmentación binaria circular — Algoritmo de segmentación

El algoritmo estándar para la segmentación.

Principio: Se toman dos puntos en la serie temporal de profundidad de los bins y se realiza una prueba t para ver si la media entre ellos es diferente. Si es significativa, se considera un punto de ruptura. Este proceso se repite.

La implementación estándar está en el paquete DNAcopy de R. CNVkit utiliza internamente este algoritmo.

Cálculo manual — Llamada de segmentos

Perfil de profundidad ficticio (razón log2).

bin 1-100: media 0,05 (referencia)
bin 101-200: media -0,95 (posible DEL heterocigota)
bin 201-300: media 0,10 (retorno a la referencia)
bin 301-500: media 0,60 (posible DUP heterocigota)
bin 501-600: media 0,00 (referencia)

Circular Binary Segmentation realizará la segmentación de la siguiente manera:

• seg 1: bin 1-100 (referencia).
• seg 2: bin 101-200 (log2=−0.95, DELE heterocigota, tamaño 100 kb).
• seg 3: bin 201-300 (referencia).
• seg 4: bin 301-500 (log2=+0.60, DUP heterocigota, tamaño 200 kb).
• seg 5: bin 501-600 (referencia).

Se determina el tipo y el tamaño de la CNV para cada segmento.

Visualización — IGV y UCSC

Al cargar el VCF de CNV en IGV, cada segmento se muestra como una caja coloreada. También es posible compararlo con la pista de número de copias del Genome Browser de UCSC.

Práctica pequeña con R en Colab

CNVkit es pesado para ejecutar en Colab; por lo tanto, es más práctico demostrar únicamente la segmentación mediante DNAcopy en R.

install.packages("DNAcopy")
library(DNAcopy)

# Datos de profundidad simulados
set.seed(42)
n <- 500
depth <- c(rnorm(100, 0, 0.1),
           rnorm(100, -1.0, 0.1),
           rnorm(100, 0, 0.1),
           rnorm(200, 0.6, 0.1))
pos <- 1:n

cna <- CNA(genomdat = depth, chrom = rep("chr1", n), maploc = pos)
cna.smoothed <- smooth.CNA(cna)
seg <- segment(cna.smoothed)

print(seg$output)
plot(seg, plot.type = "s")

Este código permite familiarizarse con la segmentación CBS y su visualización.

Mapeo de CS

• Segmentación de señales unidimensionales: consulte el episodio de DryBench «Segmentación de señales» (signal-segmentation).
• Corrección de sesgos: el contenido GC y la mapeabilidad introducen sesgos sistemáticos que pueden predecirse mediante regresión local.
• Estimación de estados basada en HMM: algunas herramientas de CNV, como gCNV, estiman con HMM los estados del número de copias (0, 1, 2, 3, 4+).

Conclusión

Las CNV son un tipo especial de SV, pero utilizan una señal completamente distinta: la profundidad de lectura. Por eso GATK Best Practices las trata en una vía separada de las SNV y los indels. En el siguiente episodio (S16) veremos cómo organizar por haplotipos las variantes llamadas por GATK e imputar con datos poblacionales conocidos las variantes no observadas en secuenciaciones pequeñas.

Para profundizar

• Talevich E. et al. (2016), CNVkit: Genome-wide copy number detection and visualization from targeted DNA sequencing. PLoS Comput Biol 12:e1004873. Artículo original de CNVkit.
• Boeva V. et al. (2012), Control-FREEC: a tool for assessing copy number and allelic content using next-generation sequencing data. Bioinformatics 28:423.
• Klambauer G. et al. (2012), cn.mops: mixture of Poissons for discovering copy number variations in next-generation sequencing data with a low false discovery rate. Nucleic Acids Research 40:e69.
• Olshen A.B. et al. (2004), Circular binary segmentation for the analysis of array-based DNA copy number data. Biostatistics 5:557. Artículo original de CBS.
• UCSC Genome Browser (https://genome.ucsc.edu/): referencia estándar para visualizar CNV.

Abrir un informe de CNVkit ayuda a comprender cómo se presenta el número de copias por gen en un informe clínico. En el siguiente episodio pasaremos al haplotipado y la imputación.