16S-V4
Para el ejercicio de caracterización de amplicones de 16S, tenemos un set de datos obtenidos de heces de rata que fueron tomadas a diferentes días post destete, algunas pocos días (early) y otras mucho después (late). Estos datos son los usados por el Dr. Patrick Schloss para su Standard Operating Procedure (SOP), publicado por Kozich et al. 2013.
El archivo a usar se llama MiSeqSOP.tar.gz (34 Mb), ya incluye un archivo de metadatos (metadata.tsv) y otro con la información del experimento (data.info), las secuencias se encuentran en una subcarpeta /seqs.
Se puede descargar directamente del servidor biobacter de la ruta /db2/curso/datasets/MiSeqSOP.tar.gz o bien de Google Drive dando click a la siguiente liga o bien desde Figshare. Para descargar desde el servidor a la carpeta donde se encuentren, se puede usar scp con el siguiente comando, cambiar user por su nombre de usuario:
$ scp user@187.141.151.196:/db2/curso/datasets/MiSeqSOP.tar.gz .
Nota. El comando termina con un punto, lo que quiere decir que descargará a la carpeta en donde están, si quieren descargarlo a otra carpeta poner la ruta, p.ej. para descargar a la carpeta Documents
$ scp user@187.141.151.196:/db2/curso/datasets/MiSeqSOP.tar.gz ~/Documents
Tabla de metadatos
La tabla de metadatos tiene un formato delimitado por tabuladores:
sample dpw time
F3D1 1 Early
F3D141 141 Late
F3D142 142 Late
...
Limpieza y ensamble
Lo primero que tenemos que hacer es limpiar y ensamblas las secuencias, para esto usaremos el pipeline del laboratorio mg_pipeline, pero antes de eso hay que descomprimirlo:
$ tar xzf MiSeqSOP.tar.gz
lo que creará una nueva carpeta llamada MiSeqSOP con una subcarpeta con las secuencias, información del estudio y una tabla de metadatos. Para limpiar las secuencias, tenemos que entrar a la carpeta seqs y ahí ejecutar el script pair-end_cleaner:
$ cd MiSeqSOP/seqs
$ pair-end_cleaner
Nos saldrá un menú con la siguiente información:
--- pair-end_cleaner v1.0.1 -------------------------
Script to unzip, clean, assemble, convert, and rename
Illumina pair-end fastq files in all 19 subdirectories.
Select a region for
16S
1 V3
2 V4
3 V3-V4
18S
4 V9
x exit
Debemos elegir la región del gen que se secuenció, para este ejemplo la opción 2, region V4 del gen 16S.
El script limpiará las secuencias para eliminar las de baja calidad, inferior a Q20, quitar las bases indefinidas n y eliminar las secuencias cortas con cutadapt. Posteriormente ensamblará las secuencias pareadas con PEAR en una solo consenso y la convertirá de fastq a fasta. Creará un archivo por para cada muestra con terminación .fna en una nueva carpeta. Esto puede tomar un par de minutos.
Eliminación de quimeras
A las secuencias limpias y ensambladas habrá que eliminarle las secuencias quiméricas generadas probablemente durante la PCR, para esto tenemos otro script basado en VSEARCH que utiliza un base de datos para este fin.
Hay que entrar a la carpeta generada en el paso anterior y que empieza con assembled.fecha_hora y ejecutar el script:
$ chimera_detector *.fna
Después de varios minutos (15 a 20), se creará una nueva subcarpeta llamada chimera_detector.fecha_hora con las secuencias libre de quimeras con terminación *.fasta
Clasificación de secuencias
Por último paso, tenemos que clasificar las secuencias para asignarles una taxonomía, esto lo haremos con el script mg_classifier:
$ mg_classifier *.fasta
Nos saldrá un menu para seleccionar la base de datos a usar
Select a database:
16S_rRNA
1 EzBioCloud-LGM 2018
2 EzBioCloud v1.5
3 SILVA v132
4 RDP 11.5
5 RDP 11.5_V3-V4
18S_rRNA
6 PR2 v4.72
7 SILVA v128_euk
x exit
En la imagen virtual anvioVM, por cuestiones de espacio, no tenemos todas las bases de datos, elijamos la número 1. Después de pocos minutos, tendremos una nueva subcarpeta con los resultados. El documento principal es una tabla de OTUs en la que está cuantas secuencias de cada OTU (líneas) se asignaron a cada muestra (columnas) analizada.
Visualización
Con R Studio podemos visualizar los resultados fácilmente siguiendo las indicaciones de este link