06.Modulo_6_ensamblaje.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Módulo 6: Ensamblajes *De novo* y control de calidad\n",
    "## Descripción general\n",
    "\n",
    "Después de generar lecturas de secuenciación y llevar a cabo una evaluación de calidad, el siguiente paso es determinar cómo encajan las lecturas buscando superposiciones entre ellas; este proceso se denomina ensamblaje del genoma. Aquí puede ver los pasos generales de un flujo de trabajo de ensamblaje del genoma.*Los pasos 1 y 2 ya se han completado.*\n",
    "\n",
    "\n",
    "![steps](images/assem.jpg)\n",
    "\n",
    "*Tomado de: https://www.ncbi.nlm.nih.gov/pmc/articles/PMC5850084/*\n",
    "\n",
    "Los datos de secuenciación por Illumina vienen como tramos relativamente cortos (35 - 150 pares de bases) de ADN, alrededor de 6 mil millones de ellos. Estas secuencias individuales se denominan lecturas de secuenciación. Hay una variedad de programas de ensamblaje que se han diseñado específicamente para ensamblar genomas a partir de datos de lectura corta de secuencias. El ensamblaje del genoma utilizando lecturas cortas (~100 pb) es complicado debido a la alta frecuencia de repeticiones más largas que la longitud de lectura de secuencia en los genomas, por ejemplo: elementos IS, operones de ARNr; y la gran cantidad de datos que los programas tienen que manejar. Además de encontrar superposiciones en la secuencia, los programas de ensamblaje también pueden usar información del tamaño de inserción predicho donde se usan lecturas emparejadas, para vincular y posicionar lecturas en un ensamblaje.\n",
    "\n",
    "Cuando un genoma se ensambla sin ninguna secuencia de referencia con la que compararlo o contra el que construir un andamiaje, se denomina ensamblaje *de novo*. El ensamblaje *de novo* puede no producir genomas completos, sino que se fragmentará en múltiples secuencias contiguas (contigs), cuyo orden es arbitrario y no guarda necesariamente ninguna relación con su orden real en el genoma.\n",
    "\n",
    "Cuando se dispone de una secuencia de referencia estrechamente relacionada, es posible mejorar un ensamblaje ordenando los contigs en comparación con la referencia y también transfiriendo la anotación. En este caso, estará presente casi todo el genoma, es decir, genes y características, pero habrá algunas regiones que contendrán lagunas o contigs que no se colocarán con precisión porque no están presentes en la referencia utilizada. Aunque técnicamente incompletos, los ensamblajes ordenados del genoma pueden proporcionar información valiosa sobre la genética y la biología de un organismo.\n",
    "\n",
    "![denovo](images/denovo.jpg)\n",
    "\n",
    "*Tomado de: https://doi.org/10.1038/nmeth.1935*\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Instalar condacolab"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!pip install -q condacolab\n",
    "import condacolab\n",
    "condacolab.install()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Instalar programas"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Instalar SPAdes, QUAST y seqtk\n",
    "!conda install bioconda::spades\n",
    "!conda install bioconda::quast\n",
    "!conda install bioconda::seqtk\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compruebe si SPAdes se instaló correctamente\n",
    "!spades.py --help"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compruebe si QUAST se instaló correctamente\n",
    "!quast.py --help"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Compruebe si seqtk se instaló correctamente\n",
    "!seqtk seq -h"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Descargar datos"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!wget https://zenodo.org/records/14231070/files/Module_6.tar.gz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Extraer el archivo .tar.gz"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!tar xvf Module_6.tar.gz"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Parte 1: Generación de un ensamblaje *de novo* para una sola cepa\n",
    "\n",
    "En esta sección, usaremos [SPAdes](https://github.com/ablab/spades), que es uno de varios ensambladores de novo que usan conjuntos de lectura cortos como entrada (por ejemplo, Illumina Reads), y el método de ensamblaje se basa en gráficos de Bruijn. \n",
    "\n",
    "Ensamblaremos el genoma de *Streptococcus pneumoniae*  \n",
    "\n",
    "Navegaremos a la carpeta que contiene los archivos pareados FASTQ ERR1795461_1.fastq.gz y ERR1795461_2.fastq.gz los cuales corresponden al run accession ERR1795461 del proyecto [PRJEB3084](https://www.ebi.ac.uk/ena/browser/view/PRJEB3084)\t\n",
    "\n",
    "Algunos datos importantes de la muestra:\n",
    "\n",
    "- País de origen: Brasil\n",
    "- Organismo: *Streptococcus pneumoniae*\n",
    "- Instrumento: ILLUMINA\n",
    "- Modelo del instrumento: Illumina MiSeq \n",
    "- Conteo de reads: 3627822\n",
    "- Conteo de bases: 453477750 \n",
    "- Nombre del centro: Wellcome Sanger Institute;SC \n",
    "- Diseño de librerías:  PAIRED\n",
    "- Estrategia de librerías: WGS"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%cd Module_6"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ahora vamos a generar una submuestra de cada archivo de secuenciación, dado que son archivos grandes y el ensamblaje completo podría requerir mucho tiempo. Sin embargo, es importante recordar que lo ideal sería realizar el análisis con la totalidad de los reads."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Submuestreo de las lecturas \n",
    "!seqtk sample -s100 ERR1795461_1.fastq.gz 100000 > ERR1795461_1_sub.fastq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Submuestreo de las lecturas\n",
    "!seqtk sample -s100 ERR1795461_2.fastq.gz 100000 > ERR1795461_2_sub.fastq"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Compresión de las lecturas\n",
    "!gzip ERR1795461_1_sub.fastq ERR1795461_2_sub.fastq"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ejecute el comando en la terminal para ejecutar SPAdes:\n",
    "\n",
    "Este paso tardará algunos minutos."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Correr SPAdes\n",
    "!spades.py -1 ERR1795461_1_sub.fastq.gz -2 ERR1795461_2_sub.fastq.gz --careful --cov-cutoff auto -o spades_assembly"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Una explicación de este comando es la siguiente:\n",
    "\n",
    "`spades.py` es la herramienta\n",
    "\n",
    "`-1` indicador para el archivo de entrada de lecturas forward\n",
    "\n",
    "`-2` indicador para el archivo de entrada de lecturas reverse\n",
    "\n",
    "`--careful` minimiza los desajustes y los indeles cortos\n",
    "\n",
    "`--cov-cutoff auto` calcula el umbral de cobertura (en lugar de la configuración predeterminada, \"apagado\")\n",
    "\n",
    "`-o` directorio de salida"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vaya al directorio de salida (spades_assembly) y mire el archivo contigs.fa:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Cambiar al directorio de salida\n",
    "%cd spades_assembly"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Consulte el informe"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Generación de un ensamblaje *de novo* para múltiples cepas\n",
    "\n",
    "También es posibles ensamblar genomas de múltiples cepas al tiempo. \n",
    "\n",
    ">**Nota:** En este módulo no ejecutaremos el ensamblaje múltiple por la falta de recursos en Colab. Sin embargo, aquí tiene un ejemplo de como hacerlo."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Primero crearemos una carpeta para cada par de archivos FASTQ comprimidos y el nombre de la identificación de la cepa usando el comando:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# No ejecutar\n",
    "# !for x in *1.fastq.gz; do mkdir ${x%%_1.fastq.gz} ; mv $x ${x%%_1.fastq.gz}; mv ${x%%1.fastq.gz}2.fastq.gz ${x%%_1.fastq.gz}; done"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Una explicación de este comando es la siguiente:\n",
    "\n",
    "`for x in *1.fastq.gz`; Esto inicia un ciclo for que itera sobre todos los archivos en el directorio actual que terminan con \"1.fastq.gz\".\n",
    "\n",
    "`do`  Esto inicia el bloque de código que se ejecutará para cada archivo.\n",
    "\n",
    "`mkdir ${x%%_1.fastq.gz}` Esto crea un nuevo directorio con el mismo nombre que el archivo, pero con el sufijo \"_1.fastq.gz\" eliminado.\n",
    "\n",
    "`mv $x ${x%%_1.fastq.gz}` Esto mueve el archivo original al nuevo directorio creado en el paso anterior.\n",
    "\n",
    "`mv \"${x%%1.fastq.gz}2.fastq.gz\" \"${x%%_1.fastq.gz}\"` Esto mueve un segundo archivo que tiene el mismo prefijo que el primer archivo, pero con un \"2.fastq sufijo .gz\", en el nuevo directorio creado en el primer paso.\n",
    "\n",
    "`done` Esto finaliza el bucle for.\n",
    "\n",
    "En general, este script está diseñado para tomar datos de secuenciación pareados que se almacenan en dos archivos separados con nombres que terminan en **\"_1.fastq.gz\"** y **\"_2.fastq.gz\"**, y organizarlos en directorios según el prefijo de los nombres de archivo."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Después se ejecuta SPAdes usando el comando:\n",
    "\n",
    ">**Nota:** Recuerda que en este módulo no ejecutaremos en ensamblaje múltiple por la falta de recursos en Colab. Sin embargo, aquí tienes un ejemplo de como hacerlo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# No ejecutar\n",
    "#!for x in *#* ; do spades.py --pe1-1 $x/${x}_1.fastq.gz --pe1-2 $x/${x}_2.fastq.gz --careful --cov-cutoff auto -o $x\"_output\"; done"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Una explicación de este comando es la siguiente:\n",
    "\n",
    "`for x in *` Esto inicia un ciclo que itera sobre todos los directorios en el directorio de trabajo actual. La variable de bucle x se establecerá en cada nombre de directorio a su vez.\n",
    "\n",
    "`do` Esta palabra clave inicia el bloque de comandos que se ejecutará para cada iteración del bucle.\n",
    "\n",
    "`spades.py` Este es el comando de ejecución SPAdes.\n",
    "\n",
    "`--pe1-1 $x/${x}_1.fastq.gz` Esto especifica la ruta al primer archivo de lectura de extremo emparejado para SPAdes. La variable $x se usa para hacer referencia al nombre del directorio actual y se agrega ${x}_1.fastq.gz para crear el nombre completo del archivo.\n",
    "\n",
    "`--pe1-2 $x/${x}_2.fastq.gz` Esto especifica la ruta al segundo archivo de lectura de extremo emparejado para SPAdes. Se usa la misma variable de nombre de directorio $x que la anterior, pero con ${x}_2.fastq.gz agregado para crear el nombre de archivo completo para el segundo archivo de lectura.\n",
    "\n",
    "`--careful` Esto le dice a SPAdes que use un algoritmo más conservador para la corrección de errores.\n",
    "\n",
    "`--cov-cutoff` auto Esto le dice a SPAdes que determine automáticamente el límite de cobertura para filtrar contigs de baja cobertura durante el ensamblaje.\n",
    "\n",
    "`-o $x\"_output\"`  Esto especifica el directorio de salida para el ensamblado de SPAdes. El nombre del directorio se crea agregando \"_output\" al nombre del directorio actual $x.\n",
    "\n",
    "`;` Este es el separador de comandos que le dice a Bash que ejecute el comando anterior antes de pasar al siguiente.\n",
    "\n",
    "`done` Esta palabra clave marca el final del bloque de bucle."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "___"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Parte 2: Control calidad para genomas ensamblados\n",
    "\n",
    "Las tecnologías modernas de secuenciación de ADN no pueden producir la secuencia completa de un cromosoma. En cambio, generan una gran cantidad de lecturas, que van desde docenas hasta miles de bases consecutivas, muestreadas de diferentes partes del genoma. Los software de ensamblaje combinan las lecturas en regiones más grandes llamadas contigs. Sin embargo, las tecnologías y softwares de secuenciación actuales enfrentan muchos desafíos que impiden la reconstrucción de cromosomas completos, incluyendo errores en las lecturas y grandes repeticiones en el genoma.\n",
    "\n",
    "Diferentes programas ensambladores utilizan diferentes enfoques heurísticos para abordar estos desafíos, lo que da como resultado muchas diferencias en los contigs que generan. Esto lleva a las preguntas de cómo evaluar la calidad de un ensamblaje y cómo comparar diferentes ensamblajes.\n",
    "\n",
    "*Otras lecturas:* https://academic.oup.com/bioinformatics/article/29/8/1072/228832"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Evaluación de ensamblajes generados\n",
    "\n",
    "Usaremos la herramienta QUAST en esta sección.\n",
    "\n",
    "[QUAST](https://github.com/ablab/quast) significa (QUality Assessment Tool) Herramienta de evaluación de calidad. Evalúa ensamblajes de genoma/metagenoma mediante el cálculo de varias métricas. El kit de herramientas QUAST actual incluye la herramienta QUAST general para ensamblajes de genomas, MetaQUAST, la extensión para conjuntos de datos metagenómicos, QUAST-LG, la extensión para genomas grandes (p. ej., mamíferos) e Icarus, el visualizador interactivo para estas herramientas.\n",
    "\n",
    "El paquete QUAST funciona con y sin genomas de referencia. Sin embargo, es mucho más informativo si se proporciona al menos un genoma de referencia cercano junto con los ensamblajes. La herramienta acepta múltiples ensamblajes, por lo que es adecuada para la comparación.\n",
    "\n",
    "Evaluaremos la calidad de los contigs generados a partir del ensamblaje *De novo* de los archivos **COMPLETOS** ERR1795461_1.fastq.gz y ERR1795461_1.fastq.qz. Es importante destacar que no se analizará la submuestra, sino un ensamblaje completo realizado previamente.\n",
    "\n",
    "Para esto, diríjase a la carpeta ERR1795461_spades_complete."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%cd .."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Ejecutaremos la herramienta QUAST en estos contigs usando el comando:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "!quast.py ERR1795461_spades_complete/contigs.fasta -r Reference_sequence_GPSC46.fa -g PROKKA_03052023.gff -1 ERR1795461_1.fastq.gz -2 ERR1795461_2.fastq.qz -o quast_ERR1795461_output"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "An explanation of this command is as follows:\n",
    "\n",
    "`quast.py` es la herramienta\n",
    "\n",
    "`-r Reference_sequence_GPSC46.fa` especifica la secuencia de referencia\n",
    "\n",
    "`-g PROKKA_03052023.gff` especifica el gen en el genoma de referencia (salida PROKKA de Reference_sequence_GPSC46.fa)\n",
    "\n",
    "`-1 ERR1795461_1.fastq.gz` archivo de entrada de lecturas forward\n",
    "\n",
    "`-2 ERR1795461_2.fastq.gz` archivo de entrada de lecturas reverse\n",
    "\n",
    "`-o quast_ERR1795461_output` especifica la carpeta de salida"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Cuando QUAST esté completo, imprimirá el siguiente mensaje en la pantalla de su terminal:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Navegue a la carpeta de salida `\"quast_ERR1795461_sub_output\"` usando el comando `%cd` y explore su contenido. Una descripción de los archivos de salida es la siguiente:\n",
    "\n",
    "`report.txt`  versión separada por tabuladores del resumen, adecuada para hojas de cálculo (Google Docs, Excel, etc.)\n",
    "\n",
    "`report.tex`  versión LaTeX del resumen\n",
    "\n",
    "`icarus.html` menú principal de Icarus con enlaces a visores interactivos.\n",
    "\n",
    "`report.pdf` todos los demás gráficos combinados con todas las tablas (el archivo se crea si está instalada la biblioteca python matplotlib)\n",
    "\n",
    "`report.html` versión HTML del informe con diagramas interactivos dentro\n",
    "\n",
    "`contigs_reports/` (solo si se proporciona un genoma de referencia)\n",
    "\n",
    "`misassemblies_report` informe detallado sobre desensamblajes\n",
    "\n",
    "`unaligned_report` informe detallado sobre contigs desalineados y parcialmente desalineados\n",
    "\n",
    "`k_mer_stats` (solo si se especifica la opción --k-mer-stats)\n",
    "\n",
    "`kmers_report`  informe detallado sobre métricas basadas en k-mer\n",
    "\n",
    "`reads_stats/` (solo si se proporcionan lecturas)\n",
    "\n",
    "`reads_report` informe detallado sobre las estadísticas de lecturas asignadas."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Interpretación\n",
    "\n",
    "### Descripción de las métricas \n",
    "\n",
    "#### 1. Informe resumido: basado en el archivo report.txt en la carpeta de salida\n",
    "\n",
    "`# contigs` es el número total de contigs en el ensamblaje.\n",
    "\n",
    "`Largest contig` es la longitud del contig más largo del ensamblaje.\n",
    "\n",
    "`Total length` es el número total de bases en el ensamblaje.\n",
    "\n",
    "`Reference length` es el número total de bases en el genoma de referencia.\n",
    "\n",
    "`GC (%)` es el número total de nucleótidos G y C en el ensamblaje, dividido por la longitud total del ensamblaje.\n",
    "\n",
    "`Reference GC (%)`  es el porcentaje de nucleótidos G y C en el genoma de referencia.\n",
    "\n",
    "`N50` es la longitud para la cual la colección de todos los contigs de esa longitud o más cubre al menos la mitad de un conjunto.\n",
    "\n",
    "`NG50` es la longitud para la cual la colección de todos los cóntigos de esa longitud o más cubre al menos la mitad del genoma de referencia (esta métrica se calcula solo si se proporciona el genoma de referencia).\n",
    "\n",
    "`Nx and NGx` (para x entre 0 y 100) ase definen de manera similar a N50 pero con x % en lugar de 50 %. El valor de x se establece con --x-for-Nx (90 por defecto).\n",
    "\n",
    "`L50 (Lx, LG50, LGx)` es el número de contigs igual o mayor que N50 (Nx, NG50, NGx)\n",
    "\n",
    "En otras palabras, L50, por ejemplo, es el número mínimo de contigs que cubren la mitad del conjunto.\n",
    "\n",
    "`# misassemblies` es el número de posiciones en los contigs (puntos de interrupción) que satisfacen uno de los siguientes criterios:\n",
    "\n",
    "- La secuencia flanqueante izquierda se alinea a más de 1 kbp de distancia de la secuencia flanqueante derecha en la referencia;\n",
    "- Las secuencias flanqueantes se superponen en más de 1 kpb;\n",
    "- Las secuencias flanqueantes se alinean con diferentes hebras o diferentes cromosomas;\n",
    "- Las secuencias flanqueantes se alinean en diferentes genomas de referencia (solo MetaQUAST).\n",
    "\n",
    "*Esta métrica requiere un genoma de referencia.*\n",
    "\n",
    "`# misassembled contigs` es el número de contigs que contienen eventos de ensamblaje incorrecto \n",
    "\n",
    "`Misassembled contigs length`  es el número total de bases en contigs mal ensamblados.\n",
    "\n",
    "`# local misassemblies` es el número de posiciones en los contigs (puntos de interrupción) que satisfacen las siguientes condiciones:\n",
    "\n",
    "1. El espacio o superposición entre las secuencias flanqueantes izquierda y derecha es inferior a 1 kpb y superior a 200 pb (la longitud máxima de indel).\n",
    "2. Las secuencias flanqueantes izquierda y derecha están ambas en la misma hebra del mismo cromosoma del genoma de referencia.\n",
    "\n",
    "`# scaffold gap ext. mis.` es el número de posiciones en los andamios (puntos de interrupción) donde las secuencias laterales se combinan en el andamio en la distancia incorrecta (suficiente para informar un gran error de montaje).\n",
    "\n",
    "`# structural variations` es el número de desensamblajes que coinciden con las variaciones estructurales del genoma (si se proporcionan lecturas o un archivo BEDPE con SV, consulte --reads1/reads2 y --sv-bedpe).\n",
    "\n",
    "`# unaligned mis. contigs`  es el número de contigs que tienen un número de bases no alineadas superior al 50 % de la longitud del contig y al menos un evento de montaje incorrecto en su fragmento alineado. Tales contigs probablemente no estén relacionados con el genoma de referencia, por lo que sus desensamblajes pueden no ser errores reales sino diferencias entre el organismo ensamblado y el de referencia.\n",
    "\n",
    "`# unaligned contigs` es el número de contigs que no tienen alineación con la secuencia de referencia. El valor \"parte X + Y\" significa X contigs totalmente desalineados más Y contigs parcialmente desalineados. Esta métrica resume # errores no alineados. contigs descritos anteriormente.\n",
    "\n",
    "`Unaligned length` es la longitud total de todas las regiones no alineadas del ensamblaje (suma de las longitudes de los contigs completamente no alineados y las partes no alineadas de las parcialmente no alineadas). Esta longitud no incluye bases no llamadas (N) en el ensamblaje.\n",
    "\n",
    "`Genome fraction (%)`  es el porcentaje de bases alineadas en el genoma de referencia. Una base en el genoma de referencia está alineada si hay al menos un cóntigo con al menos una alineación con esta base.\n",
    "\n",
    "`Duplication ratio` es el número total de bases alineadas en el ensamblaje dividido por el número total de bases alineadas en el genoma de referencia. Si el ensamblaje contiene muchos contigs que cubren las mismas regiones de la referencia, su relación de duplicación puede ser mucho mayor que 1. Esto puede ocurrir debido a la sobreestimación de multiplicidades repetidas y debido a pequeñas superposiciones entre contigs, entre otras razones.\n",
    "\n",
    "`# N's per 100 kbp` es el número promedio de bases no llamadas (N) por 100,000 bases de ensamblado.\n",
    "\n",
    "`# mismatches per 100 kbp` es el número promedio de desajustes por 100 000 bases alineadas en el ensamblaje. Los verdaderos SNP y los errores de secuenciación no se distinguen y se cuentan por igual.\n",
    "\n",
    "`# indels per 100 kbp` es el número promedio de indels por 100 000 bases alineadas en el ensamblado. Varias indels consecutivas de un solo nucleótido se cuentan como una indel.\n",
    "\n",
    "`# genomic features` es el número de características genómicas (genes, CDS, etc.) en el ensamblaje (completo y parcial), según una lista proporcionada por el usuario de posiciones de características genómicas en el genoma de referencia. Una función está \"parcialmente cubierta\" si el ensamblaje contiene al menos 100 pb de esta función, pero no la totalidad.\n",
    "\n",
    "`Total aligned length` es el número total de bases alineadas en el ensamblaje. Un valor suele ser más pequeño que un valor de longitud total porque algunos de los contigs pueden estar desalineados o parcialmente desalineados.\n",
    "\n",
    "`Largest alignment` es la longitud de la alineación continua más grande del ensamblaje. Un valor puede ser más pequeño que el valor del contig más grande si el contig más grande está mal ensamblado o parcialmente desalineado.\n",
    "\n",
    "`NA50, NGA50, NAx, NGAx, LA50, LAx, LGA50, LGAx (A\" significa \"alineado\")` son similares a las métricas correspondientes sin \"A\", pero en este caso se consideran bloques alineados en lugar de contigs. Los bloques alineados se obtienen rompiendo contigs en eventos de desmontaje y eliminando todas las bases no alineadas."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### 2. Informe resumido- descripción del gráfico (basado en icarus.html)\n",
    "\n",
    "[Icarus](https://quast.sourceforge.net/icarus) genera una visualización de tamaño de contig y una o más visualizaciones de alineación de contig (si se proporciona el genoma o los genomas de referencia). Todos ellos se encuentran en la carpeta de salida \"quast_ERR1795461_sub_output/icarus_viewers/\". Los enlaces a las visualizaciones y otra información auxiliar se proporcionan en el menú principal de Icarus, que se guarda en \"quast_ERR1795461_sub_output/icarus.html\". Todas las visualizaciones con Icarus contienen una leyenda con una descripción del esquema de colores. Para mover y hacer zoom en la ventana interactiva, puede usar el mouse, los controles de Icarus (panel superior) o los atajos de teclado (+, -, ←, →, use Shift para acelerar la acción).\n",
    "\n",
    "**Paso 1: abrir icarus.html**\n",
    "\n",
    "**Paso 2: Haga clic en contener visor de tamaño**\n",
    "\n",
    "Este tipo de visor dibuja contigs ordenados del más largo al más corto. Este orden es adecuado para comparar solo los contigs más grandes o el número de contigs más largos que un umbral específico. El visor muestra el [N50](https://quast.sourceforge.net/docs/manual.html#N50) y [Nx](https://quast.sourceforge.net/docs/manual.html#Nx) (para el valor x definido por el usuario) con color e indicación textual. Si el genoma de referencia está disponible o al menos se conoce la longitud aproximada del genoma, [NG50](https://quast.sourceforge.net/docs/manual.html#NG50) y [NGx](https://quast.sourceforge.net/docs/manual.html#Nx) también se muestran. También puede atenuar contigs más cortos que un umbral específico utilizando el panel de control de Icarus.\n",
    "\n",
    "\n",
    "**Paso 3: haz click en el reporte de QUAST en el menu principal**\n",
    "\n",
    "**Paso 4: haga clic en el visor de alineación de contig en el menú principal para ver los contig alineados con la secuencia de referencia** \n",
    "\n",
    "Este tipo de visor está disponible solo si se proporciona un genoma de referencia. Para genomas grandes (≥ 50 Mbp), cada cromosoma se muestra en un visor independiente. El espectador coloca contigs de acuerdo con su mapeo al genoma de referencia. El espectador también puede visualizar genes, operones y leer la distribución de cobertura a lo largo del genoma, si alguno de ellos se alimentó a QUAST."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Lecturas adicionales**\n",
    "\n",
    "https://quast.sourceforge.net/docs/manual.html#sec3.4 \n",
    "\n",
    "https://github.com/ablab/quast\n",
    "\n",
    "https://academic.oup.com/bioinformatics/article/29/8/1072/228832\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## BONO!\n",
    "\n",
    "Si está trabajando con BASH en su computadora o en un HPC y tiene demasiados archivos, puede optimizar los comandos, los bucles (for loops) son muy útiles para grandes conjuntos de datos.\n",
    "\n",
    "Aquí hay una manera de hacerlo.\n",
    "\n",
    "Cree un nuevo script bash usando nano llamado `assemblies.sh`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# No ejecute esta celda\n",
    "# Crea un script para ejecutar el ensamblaje\n",
    "!nano assemblies.sh"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Luego copia y pega el siguiente scrip en su nuevo archivo:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#!/bin/bash\n",
    "#Autor: Nathalia Portilla\n",
    "\n",
    "for i in $(ls *_1.trimmed.fastq.gz); do\n",
    "NAME=$(basename $i _1.trimmed.fastq.gz)\n",
    "echo \"$NAME\"\n",
    "j=\"${NAME}_1.trimmed.fastq.gz\"\n",
    "echo \"$j\"\n",
    "k=\"${NAME}_2.trimmed.fastq.gz\"\n",
    "echo \"$k\"\n",
    "spades.py -1 $j -2 $k -t 16 --careful --cov-cutoff auto -o ${NAME};\n",
    "done"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Guarda el archivo.\n",
    "\n",
    "Finalmente, puedes ejecutarlo:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# No ejecute esta celda\n",
    "!bash assemblies.sh"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Adaptado de:*\n",
    "\n",
    "- Advanced Bioinformatics Course developed for the GPS and JUNO projects - Wellcome Sanger Insitute\n",
    "    \n",
    "*Modificado por Luisa Sacristán (Universidad de los Andes-CABANA)*"
   ]
  }
 ],
 "metadata": {
  "language_info": {
   "name": "python"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}