Preparación secuencias

Antes de empezar con prácticamente cualquier análisis bioinformático con secuencias de ADN es necesario prepararlas para ello. Esta preparación puede incluir su limpieza, renombramiento de secuencias y archivos, unión de archivos o de secuencias, etc.

Limpieza de secuencias

Casi siempre el primer paso es limpiar las secuencias de adaptadores, bases indeterminadas, secuencias chicas, bases de mala calidad y otros artefactos como secuencias quiméricas.

Genómicas y metagenómicas shotgun

Si tenemos secuencias de genomas o metagenomas podemos limpiarlas dependiendo del tipo de plataforma de secuenciación aquí.

Amplicones

Sí partimos de secuencias de amplicones (16S o 18S) dependerá del tipo de plataforma:

Renombrar archivos en serie

A veces tenemos los nombres de los archivos de las muestras de una manera y queremos renombrarlos según otro código, esto lo podemos hacer si tenemos una archivo con una lista que contenga el nuevo nombre en una columna (list.txt), hay que tener cuidado que los nuevos nombres estén ordenados alfabéticamente pero de acuerdo a los nombre originales. Para asegurarnos, se puede crear un archivo excel con dos columna, en la primera el nombre viejo y en la segunda el nuevo nombre, ordenarlos por la primer columna, luego solo copiar la segunda columna y pegarla en un archivo de texto (SIN formato) y llamarla list.txt.

NOTA. Tener cuidado de revisar bien el ordenamiento de los nombres ya que excel no ordena igual que linux.

Esta nueva lista tendrá el nuevo nombre pero no la extensión requerida, por lo que si lo usamos así al final tendremos los archivos con el nuevo nombre pero sin la extensión necesaria; por ejemplo, si los nombres de los archivos fasta (*.fna), les agregaremos la extensión con sed:

$ sed -i 's/$/.fna/' list.txt

Dependiendo de cómo se escribió el listado de nuevos nombres, excel por ejemplo, es posible que el cambio de renglón no sea el que linux entiende (varía entre los sistemas operativos). Por lo tanto es buena idea primero corregir el formato con dos2unix:

$ dos2unix list.txt

Antes de cambiar todos los nombres, es recomendable hacer un respaldo de los archivos por si acaso:

$ tar czvf respaldo.tar.gz *.fna

Renombramiento:

$ for file in *.fna; do read line;  mv -v "${file}" "${line}";  done < list.txt

Una vez que concluye el renombramiento, tener cuidado de checar que se haya hecho correctamente.

Si se renombraron los archivos, es probable que el nombre de las secuencias no coincidan con los nombres de los archivos, por lo que será necesario renombrar cada una de las secuencias. Para esto tenemos un script que permite hacerlo simultáneamente para todos los nuevos archivo.

El script es fasta_header_renamer y está disponible en Github también.

Cambiar parte del nombre de archivos

La plataforma Illumina suele incluir en el nombre de los archivos una extensión con el número de muestra de la corrida, de la forma _S1, _S10, etc. Obviamente esta extensión nos complica luego el análisis y es conveniente eliminarla.

Por ejemplo, a la muestra llamada 1T se le añadió _S87:

├── 1T_S87_L001_R1_001.fastq
├── 1T_S87_L001_R2_001.fastq

Tenemos que eliminar _S87 de ambos archivos, pero esto se complica cuando tenemos muchas muestras y éstas tienen una extensión diferente; al menos todas tienen un formato similar (_S...), el problema es que el número consecutivo que tienen puede ser de uno o dos dígitos, y esto lo complica.

Lo más fácil es eliminar esta extensión en los archivos fasta ya limpios.

El script pair-end_cleaner ya se encarga de eliminar esta extensión automáticamente.

Manualmente podemos hacerlo con el siguiente comando:

$ rename 's/_S\d*//' *.fastq

Unión de archivos

Archivos fasta

Es común que la misma muestra se secuencíe más de una vez por lo que será necesario juntar ambos archivos antes o después de limpiarlos; por ejemplo en el siguiente directorio existen dos subdirectorios con un archivo con el mismo nombre (M1.fasta) y otros archivos con nombres diferentes (M2.fasta y M3.fasta). En este caso, la muestra M1 se secuenció dos veces, quizá porque en la primera secuenciación se obtuvieron pocas secuencias. 

.
├── folderA
│   ├── M1.fasta
│   └── M2.fasta
└── folderB
    ├── M1.fasta
    └── M3.fasta

Para este caso, es sencillo unirlos con el comando cat :

$ cat folderA/M1.fasta folderB/M1.fasta > M1_merged.fasta

Así crearemos un nuevo archivo con la unión de ambos archivos llamado M1_merged.fasta.

Pero resulta mas latoso si son muchas muestras las que es necesario concatenar y los archivos se encuentren en varios directorios. Para esto podemos usar un script que lo facilite, pero los archivos tienen por fuerza que tener el mismo nombre, lo que no sucede si se resecuenció la muestra pues en cada corrida, Illumina le añade al archivo unos caracteres únicos de la forma _S00_ así, tendríamos primero que renombrarlos como se explica arriba para eliminar esto del nombre.

En Biobacter tenemos un pequeño script (files_merger.sh) que busca archivos con el mismo nombre en subdirectorios y los concatena, además también copia los archivos que no tienen otros en un nuevo directorio junto con los concatenados. Los archivos originales no los toca.

Script para buscar archivos fasta con el mismo nombre y concatenarlos:

file_merger.sh

Archivos fastq

La unión de archivos podría hacerse incluso antes de la limpieza, usando los archivos generados directamente de las corridas Illumina pair-end. La estructura de archivos de este tipo de secuenciaciones son similares a esta:

.
├── FASTQ_Generation_2020-08-23_16_31_53Z-302060759
│   ├── SAMPLE1_L001-ds.f7d621caecda4d019180882537cb1dc6
│   │   ├── SAMPLE1_S34_L001_R1_001.fastq.gz
│   │   └── SAMPLE1_S34_L001_R2_001.fastq.gz
│   ├── SAMPLE2_L001-ds.dc55815d3139446ba4cbcbaabbc0a040
│   │   ├── SAMPLE2_S35_L001_R1_001.fastq.gz
│   │   └── SAMPLE2_S35_L001_R2_001.fastq.gz
│   └── SAMPLE3_L001-ds.d77071d7bfd640c68af16219ea8508c0
│       ├── SAMPLE3_S36_L001_R1_001.fastq.gz
│       └── SAMPLE3_S36_L001_R2_001.fastq.gz
└── FASTQ_Generation_2020-10-07_19_47_56Z-325516191
    └── SAMPLE3_L001-ds.492a4bffaed24ca9be4a6b8efc2fbdd4
        ├── SAMPLE3_S19_L001_R1_001.fastq.gz
        └── SAMPLE3_S19_L001_R2_001.fastq.gz

Hubo dos corridas de secuenciación (directorios FASTQ_Generation_ con diferente fecha), la primer corrida con tres muestras (SAMPLE1 a SAMPLE3) y en la segunda corrida se resecuenció la muestra SAMPLE3. La idea es sólo concatenar los archivos R1 y R2 de ambas corridas para ésta muestra solamente y los demás dejarlos como están, pues esos ya tienen suficientes secuencias para hacer un buen análisis.

Podríamos hacer esto fácilmente con cat, como se explicó arriba, pero cuando son muchas muestras a juntar, es muy engorroso, para esto tenemos también un script que lo hace automáticamente fastq_file_merger. Este script tiene que tener la estructura de archivos mencionada anteriormente obligatoriamente.

Script para buscar archivos fastq con el mismo nombre en diferentes directorios y concatenarlos:

fastq_file_merger.sh

Eliminar secuencias

A veces para algunos procesos necesitamos eliminar ciertas secuencias con bases ambiguas, como N, Y, K, etc.

Primero podemos ver si tenemos bases ambiguas en las secuencias en archivos multifasta (p.ej. llamado file.fasta) sin considerar el header de las mismas secuencias:

$ grep -v ">" file.fasta | grep -E -i -c 'R|Y|K|M|S|W|B|D|H|V|N'

Con el primer comando de este pipe (grep -v ">") no consideramos el header, que siempre comienza con >. Con el segundo comando buscamos las bases ambiguas solo en las secuencias. Si solo queremos ver si tenemos un solo tipo de ambigüedad, quitaríamos las otras: grep -E -i -c 'N' en este caso solo para buscar Ns.

Para eliminar solo las secuencia, incluyendo su header que contengas N:

$ awk '/^>/{print s? s"\n"$0:$0;s="";next}{s=s sprintf("%s",$0)}END{if(s)print s}' file.fasta | awk -v RS='>' '!/N/{printf $0RT}'