Saltar a contenido

Parte 3: Llamado de variantes conjunto en una cohorte

Traducción asistida por IA - más información y sugerencias

Anteriormente, construiste un pipeline de llamado de variantes por muestra que procesaba los datos de cada muestra de forma independiente. Ahora vamos a extenderlo para implementar el llamado de variantes conjunto, como se describe en la Parte 1 (caso de uso multi-muestra).

Cómo comenzar desde esta sección

Esta sección del curso asume que has completado la Parte 1: Descripción del método, la Parte 2: Llamado de variantes por muestra y tienes un pipeline genomics.nf funcional.

Si no completaste la Parte 2 o quieres comenzar de nuevo para esta parte, puedes usar la solución de la Parte 2 como punto de partida. Ejecuta estos comandos desde dentro del directorio nf4-science/genomics/:

cp solutions/part2/genomics-2.nf genomics.nf
cp solutions/part2/nextflow.config .
cp solutions/part2/modules/* modules/

Esto te da un workflow completo de llamado de variantes por muestra. Puedes probar que se ejecuta exitosamente ejecutando el siguiente comando:

nextflow run genomics.nf -profile test

Asignación

En esta parte del curso, vamos a extender el workflow para hacer lo siguiente:

  1. Generar un archivo de índice para cada archivo BAM de entrada usando Samtools
  2. Ejecutar GATK HaplotypeCaller en cada archivo BAM de entrada para generar un GVCF de llamados de variantes genómicas por muestra
  3. Recolectar todos los GVCFs y combinarlos en un almacén de datos GenomicsDB
  4. Ejecutar genotipado conjunto en el almacén de datos GVCF combinado para producir un VCF a nivel de cohorte
BAMSamtools indexBAM indexIntervalsReference+ index & dictGATK HaplotypeCallerGVCF + indexx multiple samplesGVCF + indexGVCF + indexGVCF + indexGenomicsDBvariant storeGATK GenomicsDBImportGATK GenotypeGVCFsJoint-calledVCFGVCF mode

Esto implementa el método descrito en la Parte 1: Descripción del método (segunda sección que cubre el caso de uso multi-muestra) y se construye directamente sobre el workflow producido en la Parte 2: Llamado de variantes por muestra.

Plan de la lección

Hemos dividido esto en dos etapas:

  1. Modificar el paso de llamado de variantes por muestra para producir un GVCF. Esto cubre la actualización de comandos y salidas del proceso.
  2. Agregar un paso de genotipado conjunto que combine y genotipe los GVCFs por muestra. Esto introduce el operador collect(), closures de Groovy para construcción de línea de comandos, y procesos con múltiples comandos.

Esto automatiza los pasos de la segunda sección de la Parte 1: Descripción del método, donde ejecutaste estos comandos manualmente en sus contenedores.

Consejo

Asegúrate de estar en el directorio de trabajo correcto: cd /workspaces/training/nf4-science/genomics

1. Modificar el paso de llamado de variantes por muestra para producir un GVCF

El pipeline de la Parte 2 produce archivos VCF, pero el llamado conjunto requiere archivos GVCF. Necesitamos activar el modo de llamado de variantes GVCF y actualizar la extensión del archivo de salida.

Recuerda el comando de llamado de variantes GVCF de la Parte 1:

gatk HaplotypeCaller \
        -R /data/ref/ref.fasta \
        -I /data/bam/reads_mother.bam \
        -O reads_mother.g.vcf \
        -L /data/ref/intervals.bed \
        -ERC GVCF

Comparado con el comando base de HaplotypeCaller que envolvimos en la Parte 2, las diferencias son el parámetro -ERC GVCF y la extensión de salida .g.vcf.

1.1. Indicar a HaplotypeCaller que emita un GVCF y actualizar la extensión de salida

Abre el archivo de módulo modules/gatk_haplotypecaller.nf para hacer dos cambios:

  • Agregar el parámetro -ERC GVCF al comando GATK HaplotypeCaller;
  • Actualizar la ruta del archivo de salida para usar la extensión .g.vcf correspondiente, según la convención de GATK.

Asegúrate de agregar una barra invertida (\) al final de la línea anterior cuando agregues -ERC GVCF.

modules/gatk_haplotypecaller.nf
22
23
24
25
26
27
28
29
    """
    gatk HaplotypeCaller \
        -R ${ref_fasta} \
        -I ${input_bam} \
        -O ${input_bam}.g.vcf \
        -L ${interval_list} \
        -ERC GVCF
    """
modules/gatk_haplotypecaller.nf
22
23
24
25
26
27
28
    """
    gatk HaplotypeCaller \
        -R ${ref_fasta} \
        -I ${input_bam} \
        -O ${input_bam}.vcf \
        -L ${interval_list}
    """

También necesitamos actualizar el bloque de salida para que coincida con la nueva extensión de archivo. Dado que cambiamos la salida del comando de .vcf a .g.vcf, el bloque output: del proceso debe reflejar el mismo cambio.

1.2. Actualizar la extensión del archivo de salida en el bloque de salidas del proceso

modules/gatk_haplotypecaller.nf
17
18
19
    output:
    path "${input_bam}.g.vcf"     , emit: vcf
    path "${input_bam}.g.vcf.idx" , emit: idx
modules/gatk_haplotypecaller.nf
17
18
19
    output:
    path "${input_bam}.vcf"     , emit: vcf
    path "${input_bam}.vcf.idx" , emit: idx

También necesitamos actualizar la configuración de publicación y salida del workflow para reflejar las nuevas salidas GVCF.

1.3. Actualizar los destinos de publicación para las nuevas salidas GVCF

Dado que ahora estamos produciendo GVCFs en lugar de VCFs, deberíamos actualizar la sección publish: del workflow para usar nombres más descriptivos. También organizaremos los archivos GVCF en su propio subdirectorio para mayor claridad.

genomics.nf
47
48
49
50
    publish:
    indexed_bam = SAMTOOLS_INDEX.out
    gvcf = GATK_HAPLOTYPECALLER.out.vcf
    gvcf_idx = GATK_HAPLOTYPECALLER.out.idx
genomics.nf
47
48
49
50
    publish:
    indexed_bam = SAMTOOLS_INDEX.out
    vcf = GATK_HAPLOTYPECALLER.out.vcf
    vcf_idx = GATK_HAPLOTYPECALLER.out.idx

Ahora actualiza el bloque de salida para que coincida.

1.4. Actualizar el bloque de salida para la nueva estructura de directorios

También necesitamos actualizar el bloque output para colocar los archivos GVCF en un subdirectorio gvcf.

genomics.nf
53
54
55
56
57
58
59
60
61
62
63
output {
    indexed_bam {
        path 'indexed_bam'
    }
    gvcf {
        path 'gvcf'
    }
    gvcf_idx {
        path 'gvcf'
    }
}
genomics.nf
53
54
55
56
57
58
59
60
61
62
63
output {
    indexed_bam {
        path 'bam'
    }
    vcf {
        path 'vcf'
    }
    vcf_idx {
        path 'vcf'
    }
}

Con el módulo, los destinos de publicación y el bloque de salida todos actualizados, podemos probar los cambios.

1.5. Ejecutar el pipeline

Ejecuta el workflow para verificar que los cambios funcionen.

nextflow run genomics.nf -profile test
Salida del comando 
N E X T F L O W   ~  version 25.10.2

┃ Launching `genomics.nf` [nostalgic_franklin] DSL2 - revision: f2c0a93c6a

executor >  local (6)
[cc/fbc705] SAMTOOLS_INDEX (3)       | 3 of 3 ✔
[27/0d7eb9] GATK_HAPLOTYPECALLER (2) | 3 of 3 ✔

La salida de Nextflow se ve igual que antes, pero los archivos .g.vcf y sus archivos de índice ahora están organizados en subdirectorios.

Contenido del directorio (enlaces simbólicos acortados) 
results/
├── gvcf/
│   ├── reads_father.bam.g.vcf -> */27/0d7eb9*/reads_father.bam.g.vcf
│   ├── reads_father.bam.g.vcf.idx -> */27/0d7eb9*/reads_father.bam.g.vcf.idx
│   ├── reads_mother.bam.g.vcf -> */e4/4ed55e*/reads_mother.bam.g.vcf
│   ├── reads_mother.bam.g.vcf.idx -> */e4/4ed55e*/reads_mother.bam.g.vcf.idx
│   ├── reads_son.bam.g.vcf -> */08/e95962*/reads_son.bam.g.vcf
│   └── reads_son.bam.g.vcf.idx -> */08/e95962*/reads_son.bam.g.vcf.idx
└── indexed_bam/
    ├── reads_father.bam -> */9a/c7a873*/reads_father.bam
    ├── reads_father.bam.bai -> */9a/c7a873*/reads_father.bam.bai
    ├── reads_mother.bam -> */f1/8d8486*/reads_mother.bam
    ├── reads_mother.bam.bai -> */f1/8d8486*/reads_mother.bam.bai
    ├── reads_son.bam -> */cc/fbc705*/reads_son.bam
    └── reads_son.bam.bai -> */cc/fbc705*/reads_son.bam.bai

Si abres uno de los archivos GVCF y lo revisas, puedes verificar que GATK HaplotypeCaller produjo archivos GVCF como se solicitó.

Conclusión

Cuando cambias el nombre del archivo de salida de un comando de herramienta, el bloque output: del proceso y la configuración de publicación/salida deben actualizarse para coincidir.

¿Qué sigue?

Aprende a recolectar el contenido de un canal y pasarlo al siguiente proceso como una sola entrada.

2. Agregar un paso de genotipado conjunto

Ahora necesitamos recolectar los GVCFs por muestra, combinarlos en un almacén de datos GenomicsDB y ejecutar el genotipado conjunto para producir un VCF a nivel de cohorte. Como se cubrió en la Parte 1, esta es una operación de dos herramientas: GenomicsDBImport combina los GVCFs, luego GenotypeGVCFs produce los llamados de variantes finales. Envolveremos ambas herramientas en un solo proceso llamado GATK_JOINTGENOTYPING.

Recuerda los dos comandos de la Parte 1:

gatk GenomicsDBImport \
    -V reads_mother.g.vcf \
    -V reads_father.g.vcf \
    -V reads_son.g.vcf \
    -L /data/ref/intervals.bed \
    --genomicsdb-workspace-path family_trio_gdb

gatk GenotypeGVCFs \
    -R /data/ref/ref.fasta \
    -V gendb://family_trio_gdb \
    -O family_trio.vcf

El primer comando toma los GVCFs por muestra y un archivo de intervalos, y produce un almacén de datos GenomicsDB. El segundo toma ese almacén de datos, un genoma de referencia, y produce el VCF final a nivel de cohorte. El URI del contenedor es el mismo que para HaplotypeCaller: community.wave.seqera.io/library/gatk4:4.5.0.0--730ee8817e436867.

2.1. Configurar las entradas

El proceso de genotipado conjunto necesita dos tipos de entradas que aún no tenemos: un nombre de cohorte arbitrario, y las salidas GVCF recolectadas de todas las muestras agrupadas juntas.

2.1.1. Agregar un parámetro cohort_name

Necesitamos proporcionar un nombre arbitrario para la cohorte. Más adelante en la serie de capacitación aprenderás cómo usar metadatos de muestras para este tipo de cosas, pero por ahora simplemente declaramos un parámetro CLI usando params.

genomics.nf
18
19
20
21
22
    intervals: Path

    // Nombre base para el archivo de salida final
    cohort_name: String
}
genomics.nf
18
19
    intervals: Path
}

Usaremos este parámetro para nombrar el archivo de salida final.

2.1.2. Agregar un valor predeterminado para cohort_name en el perfil de prueba

También agregamos un valor predeterminado para el parámetro cohort_name en el perfil de prueba:

nextflow.config
 4
 5
 6
 7
 8
 9
10
11
test {
    params.input = "${projectDir}/data/samplesheet.csv"
    params.reference = "${projectDir}/data/ref/ref.fasta"
    params.reference_index = "${projectDir}/data/ref/ref.fasta.fai"
    params.reference_dict = "${projectDir}/data/ref/ref.dict"
    params.intervals = "${projectDir}/data/ref/intervals.bed"
    params.cohort_name = "family_trio"
}
nextflow.config
 4
 5
 6
 7
 8
 9
10
test {
    params.input = "${projectDir}/data/samplesheet.csv"
    params.reference = "${projectDir}/data/ref/ref.fasta"
    params.reference_index = "${projectDir}/data/ref/ref.fasta.fai"
    params.reference_dict = "${projectDir}/data/ref/ref.dict"
    params.intervals = "${projectDir}/data/ref/intervals.bed"
}

A continuación, necesitaremos reunir las salidas por muestra para que puedan procesarse juntas.

2.1.3. Reunir las salidas de HaplotypeCaller a través de las muestras

Si conectáramos el canal de salida de GATK_HAPLOTYPECALLER directamente al nuevo proceso, Nextflow llamaría al proceso en cada GVCF de muestra por separado. Queremos agrupar los tres GVCFs (y sus archivos de índice) para que Nextflow los entregue todos juntos a una sola llamada de proceso.

Podemos hacer eso usando el operador de canal collect(). Agrega las siguientes líneas al cuerpo del workflow, justo después de la llamada a GATK_HAPLOTYPECALLER:

genomics.nf
48
49
50
51
52
53
        intervals_file
    )

    // Recolectar salidas de llamado de variantes a través de las muestras
    all_gvcfs_ch = GATK_HAPLOTYPECALLER.out.vcf.collect()
    all_idxs_ch = GATK_HAPLOTYPECALLER.out.idx.collect()
genomics.nf
48
49
        intervals_file
    )

Desglosando esto:

  1. Tomamos el canal de salida de GATK_HAPLOTYPECALLER usando la propiedad .out.
  2. Debido a que nombramos las salidas usando emit: en la sección 1, podemos seleccionar los GVCFs con .vcf y los archivos de índice con .idx. Sin salidas nombradas, tendríamos que usar .out[0] y .out[1].
  3. El operador collect() agrupa todos los archivos en un solo elemento, por lo que all_gvcfs_ch contiene los tres GVCFs juntos, y all_idxs_ch contiene los tres archivos de índice juntos.

Podemos recolectar los GVCFs y sus archivos de índice por separado (en lugar de mantenerlos juntos en tuplas) porque Nextflow preparará todos los archivos de entrada juntos para la ejecución, por lo que los archivos de índice estarán presentes junto a los GVCFs.

Consejo

Puedes usar el operador view() para inspeccionar el contenido de los canales antes y después de aplicar operadores de canal.

2.2. Escribir el proceso de genotipado conjunto y llamarlo en el workflow

Siguiendo el mismo patrón que usamos en la Parte 2, escribiremos la definición del proceso en un archivo de módulo, lo importaremos al workflow y lo llamaremos con las entradas que acabamos de preparar.

2.2.1. Construir una cadena para dar a cada GVCF un argumento -V

Antes de comenzar a completar la definición del proceso, hay una cosa que resolver. El comando GenomicsDBImport espera un argumento -V separado para cada archivo GVCF, así:

gatk GenomicsDBImport \
    -V reads_mother.bam.g.vcf \
    -V reads_father.bam.g.vcf \
    -V reads_son.bam.g.vcf \
    ...

Si escribiéramos -V ${all_gvcfs_ch}, Nextflow simplemente concatenaría los nombres de archivo y esa parte del comando se vería así:

-V reads_mother.bam.g.vcf reads_father.bam.g.vcf reads_son.bam.g.vcf

Pero necesitamos que la cadena se vea así:

-V reads_mother.bam.g.vcf -V reads_father.bam.g.vcf -V reads_son.bam.g.vcf

Importante, necesitamos construir esta cadena dinámicamente a partir de cualquier archivo que esté en el canal recolectado. Nextflow (vía Groovy) proporciona una forma concisa de hacer esto:

def gvcfs_line = all_gvcfs.collect { gvcf -> "-V ${gvcf}" }.join(' ')

Desglosando esto:

  1. all_gvcfs.collect { gvcf -> "-V ${gvcf}" } itera sobre cada ruta de archivo y antepone -V a ella, produciendo ["-V A.g.vcf", "-V B.g.vcf", "-V C.g.vcf"].
  2. .join(' ') las concatena con espacios: "-V A.g.vcf -V B.g.vcf -V C.g.vcf".
  3. El resultado se asigna a una variable local gvcfs_line (definida con def), que podemos interpolar en la plantilla de comando.

Esta línea va dentro del bloque script: del proceso, antes de la plantilla de comando. Puedes colocar código Groovy arbitrario entre script: y las """ de apertura de la plantilla de comando.

Luego podrás referirte a toda esa cadena como gvcfs_line en el bloque script: del proceso.

2.2.2. Completar el módulo para el proceso de genotipado conjunto

A continuación, podemos abordar la escritura del proceso completo.

Abre modules/gatk_jointgenotyping.nf y examina el esquema de la definición del proceso.

Adelante, completa la definición del proceso usando la información proporcionada arriba, luego verifica tu trabajo contra la solución en la pestaña "Después" a continuación.

modules/gatk_jointgenotyping.nf
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
#!/usr/bin/env nextflow

/*
 * Combinar GVCFs en almacén de datos GenomicsDB y ejecutar genotipado conjunto para producir llamados a nivel de cohorte
 */
process GATK_JOINTGENOTYPING {

    container

    input:

    output:

    script:
    """

    """
}
modules/gatk_jointgenotyping.nf
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
#!/usr/bin/env nextflow

/*
 * Combinar GVCFs en almacén de datos GenomicsDB y ejecutar genotipado conjunto para producir llamados a nivel de cohorte
 */
process GATK_JOINTGENOTYPING {

    container "community.wave.seqera.io/library/gatk4:4.5.0.0--730ee8817e436867"

    input:
    path all_gvcfs
    path all_idxs
    path interval_list
    val cohort_name
    path ref_fasta
    path ref_index
    path ref_dict

    output:
    path "${cohort_name}.joint.vcf"     , emit: vcf
    path "${cohort_name}.joint.vcf.idx" , emit: idx

    script:
    def gvcfs_line = all_gvcfs.collect { gvcf -> "-V ${gvcf}" }.join(' ')
    """
    gatk GenomicsDBImport \
        ${gvcfs_line} \
        -L ${interval_list} \
        --genomicsdb-workspace-path ${cohort_name}_gdb

    gatk GenotypeGVCFs \
        -R ${ref_fasta} \
        -V gendb://${cohort_name}_gdb \
        -L ${interval_list} \
        -O ${cohort_name}.joint.vcf
    """
}

Hay varias cosas que vale la pena destacar aquí.

Como anteriormente, varias entradas están listadas aunque los comandos no las referencian directamente: all_idxs, ref_index y ref_dict. Listarlas asegura que Nextflow prepare estos archivos en el directorio de trabajo junto a los archivos que sí aparecen en los comandos, que GATK espera encontrar basándose en convenciones de nomenclatura.

La variable gvcfs_line usa el closure de Groovy descrito arriba para construir los argumentos -V para GenomicsDBImport.

Este proceso ejecuta dos comandos en serie, tal como lo harías en la terminal. GenomicsDBImport combina los GVCFs por muestra en un almacén de datos, luego GenotypeGVCFs lee ese almacén de datos y produce el VCF final a nivel de cohorte. El almacén de datos GenomicsDB (${cohort_name}_gdb) es un artefacto intermedio usado solo dentro del proceso; no aparece en el bloque de salida.

Una vez que hayas completado esto, el proceso está listo para usar. Para usarlo en el workflow, necesitarás importar el módulo y agregar una llamada de proceso.

2.2.3. Importar el módulo

Agrega la declaración de importación a genomics.nf, debajo de las declaraciones de importación existentes:

genomics.nf
3
4
5
6
// Declaraciones INCLUDE de módulos
include { SAMTOOLS_INDEX } from './modules/samtools_index.nf'
include { GATK_HAPLOTYPECALLER } from './modules/gatk_haplotypecaller.nf'
include { GATK_JOINTGENOTYPING } from './modules/gatk_jointgenotyping.nf'
genomics.nf
3
4
5
// Declaraciones INCLUDE de módulos
include { SAMTOOLS_INDEX } from './modules/samtools_index.nf'
include { GATK_HAPLOTYPECALLER } from './modules/gatk_haplotypecaller.nf'

El proceso ahora está disponible en el ámbito del workflow.

2.2.4. Agregar la llamada al proceso

Agrega la llamada a GATK_JOINTGENOTYPING en el cuerpo del workflow, después de las líneas collect():

genomics.nf
53
54
55
56
57
58
59
60
61
62
63
64
    all_idxs_ch = GATK_HAPLOTYPECALLER.out.idx.collect()

    // Combinar GVCFs en un almacén de datos GenomicsDB y aplicar genotipado conjunto
    GATK_JOINTGENOTYPING(
        all_gvcfs_ch,
        all_idxs_ch,
        intervals_file,
        params.cohort_name,
        ref_file,
        ref_index_file,
        ref_dict_file
    )
genomics.nf
53
    all_idxs_ch = GATK_HAPLOTYPECALLER.out.idx.collect()

El proceso ahora está completamente conectado. A continuación, configuramos cómo se publican las salidas.

2.3. Configurar el manejo de salidas

Necesitamos publicar las salidas del VCF conjunto. Agrega destinos de publicación y entradas de bloque de salida para los resultados del genotipado conjunto.

2.3.1. Agregar destinos de publicación para el VCF conjunto

Agrega el VCF conjunto y su índice a la sección publish: del workflow:

genomics.nf
66
67
68
69
70
71
    publish:
    indexed_bam = SAMTOOLS_INDEX.out
    gvcf = GATK_HAPLOTYPECALLER.out.vcf
    gvcf_idx = GATK_HAPLOTYPECALLER.out.idx
    joint_vcf = GATK_JOINTGENOTYPING.out.vcf
    joint_vcf_idx = GATK_JOINTGENOTYPING.out.idx
genomics.nf
66
67
68
69
    publish:
    indexed_bam = SAMTOOLS_INDEX.out
    gvcf = GATK_HAPLOTYPECALLER.out.vcf
    gvcf_idx = GATK_HAPLOTYPECALLER.out.idx

Ahora actualiza el bloque de salida para que coincida.

2.3.2. Agregar entradas de bloque de salida para el VCF conjunto

Agrega entradas para los archivos VCF conjuntos. Los colocaremos en la raíz del directorio de resultados ya que esta es la salida final.

genomics.nf
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
output {
    indexed_bam {
        path 'indexed_bam'
    }
    gvcf {
        path 'gvcf'
    }
    gvcf_idx {
        path 'gvcf'
    }
    joint_vcf {
        path '.'
    }
    joint_vcf_idx {
        path '.'
    }
}
genomics.nf
74
75
76
77
78
79
80
81
82
83
84
output {
    indexed_bam {
        path 'indexed_bam'
    }
    gvcf {
        path 'gvcf'
    }
    gvcf_idx {
        path 'gvcf'
    }
}

Con el proceso, los destinos de publicación y el bloque de salida todos en su lugar, podemos probar el workflow completo.

2.4. Ejecutar el workflow

Ejecuta el workflow para verificar que todo funcione.

nextflow run genomics.nf -profile test -resume
Salida del comando 
N E X T F L O W   ~  version 25.10.2

┃ Launching `genomics.nf` [crazy_marconi] DSL2 - revision: 5da9afc841

executor >  local (1)
[9a/c7a873] SAMTOOLS_INDEX (2)       | 3 of 3, cached: 3 ✔
[e4/4ed55e] GATK_HAPLOTYPECALLER (2) | 3 of 3, cached: 3 ✔
[a6/7cc8ed] GATK_JOINTGENOTYPING     | 1 of 1 ✔

Los primeros dos pasos están en caché de la ejecución anterior, y el nuevo paso GATK_JOINTGENOTYPING se ejecuta una vez en las entradas recolectadas de las tres muestras. El archivo de salida final, family_trio.joint.vcf (y su índice), están en el directorio de resultados.

Contenido del directorio (enlaces simbólicos acortados) 
results/
├── family_trio.joint.vcf -> */a6/7cc8ed*/family_trio.joint.vcf
├── family_trio.joint.vcf.idx -> */a6/7cc8ed*/family_trio.joint.vcf.idx
├── gvcf/
│   ├── reads_father.bam.g.vcf -> */27/0d7eb9*/reads_father.bam.g.vcf
│   ├── reads_father.bam.g.vcf.idx -> */27/0d7eb9*/reads_father.bam.g.vcf.idx
│   ├── reads_mother.bam.g.vcf -> */e4/4ed55e*/reads_mother.bam.g.vcf
│   ├── reads_mother.bam.g.vcf.idx -> */e4/4ed55e*/reads_mother.bam.g.vcf.idx
│   ├── reads_son.bam.g.vcf -> */08/e95962*/reads_son.bam.g.vcf
│   └── reads_son.bam.g.vcf.idx -> */08/e95962*/reads_son.bam.g.vcf.idx
└── indexed_bam/
    ├── reads_father.bam -> */9a/c7a873*/reads_father.bam
    ├── reads_father.bam.bai -> */9a/c7a873*/reads_father.bam.bai
    ├── reads_mother.bam -> */f1/8d8486*/reads_mother.bam
    ├── reads_mother.bam.bai -> */f1/8d8486*/reads_mother.bam.bai
    ├── reads_son.bam -> */cc/fbc705*/reads_son.bam
    └── reads_son.bam.bai -> */cc/fbc705*/reads_son.bam.bai

Si abres el archivo VCF conjunto, puedes verificar que el workflow produjo los llamados de variantes esperados.

family_trio.joint.vcf
40
41
42
43
#CHROM	POS	ID	REF	ALT	QUAL	FILTER	INFO	FORMAT	reads_father	reads_mother	reads_son
20_10037292_10066351	3480	.	C	CT	1625.89	.	AC=5;AF=0.833;AN=6;BaseQRankSum=0.220;DP=85;ExcessHet=0.0000;FS=2.476;MLEAC=5;MLEAF=0.833;MQ=60.00;MQRankSum=0.00;QD=21.68;ReadPosRankSum=-1.147e+00;SOR=0.487	GT:AD:DP:GQ:PL	0/1:15,16:31:99:367,0,375	1/1:0,18:18:54:517,54,0	1/1:0,26:26:78:756,78,0
20_10037292_10066351	3520	.	AT	A	1678.89	.	AC=5;AF=0.833;AN=6;BaseQRankSum=1.03;DP=80;ExcessHet=0.0000;FS=2.290;MLEAC=5;MLEAF=0.833;MQ=60.00;MQRankSum=0.00;QD=22.39;ReadPosRankSum=0.701;SOR=0.730	GT:AD:DP:GQ:PL	0/1:18,13:31:99:296,0,424	1/1:0,18:18:54:623,54,0	1/1:0,26:26:78:774,78,0
20_10037292_10066351	3529	.	T	A	154.29	.	AC=1;AF=0.167;AN=6;BaseQRankSum=-5.440e-01;DP=104;ExcessHet=0.0000;FS=1.871;MLEAC=1;MLEAF=0.167;MQ=60.00;MQRankSum=0.00;QD=7.71;ReadPosRankSum=-1.158e+00;SOR=1.034	GT:AD:DP:GQ:PL	0/0:44,0:44:99:0,112,1347	0/1:12,8:20:99:163,0,328	0/0:39,0:39:99:0,105,1194

¡Ahora tienes un workflow de llamado de variantes conjunto automatizado y completamente reproducible!

Nota

Ten en cuenta que los archivos de datos que te dimos cubren solo una pequeña porción del cromosoma 20. El tamaño real de un conjunto de llamados de variantes se contaría en millones de variantes. ¡Por eso usamos solo subconjuntos pequeños de datos para propósitos de capacitación!

Conclusión

Sabes cómo recolectar salidas de un canal y agruparlas como una sola entrada a otro proceso. También sabes cómo construir una línea de comandos usando closures de Groovy, y cómo ejecutar múltiples comandos en un solo proceso.

¿Qué sigue?

¡Date una gran palmada en la espalda! Has completado el curso Nextflow para Genómica.

Dirígete al resumen del curso final para revisar lo que aprendiste y descubrir qué viene después.