Manejo de datos de secuenciación con Bash

Al finalizar este tema

En la terminal, podrás abrir archivos FASTA/FASTQ, contar el número de secuencias y buscar secuencias específicas. Esta es tu primera experiencia manipulando datos solo con el teclado, sin usar el ratón.

¿Qué es una terminal?

Así como en un laboratorio se manipulan con precisión los reactivos usando una pipeta, la terminal es una herramienta para manejar con precisión los datos de la computadora. En lugar de hacer clic, se introducen comandos.

Al agregar el prefijo ! a una celda en Google Colab, podrás ejecutar comandos de Bash.

Crear y examinar un archivo

Primero, crearemos un archivo FASTA para fines prácticos.

# Crear un archivo FASTA para la práctica
echo ">BRCA1_human
ATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC
>TP53_human
ATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA
>EGFR_human
ATGCGACCCTCCGGGACGGCCGGGGCAGCGCTCCTGGCGCTG" > genes.fasta

# Comprobar el contenido del archivo
cat genes.fasta

# Contar las líneas del archivo
wc -l genes.fasta

wc es la abreviatura de "Word Count". La opción -l cuenta solo el número de líneas.

grep: Buscar secuencias

grep es un comando que busca patrones específicos en un archivo, de forma similar a la búsqueda de palabras clave en un artículo científico.

# Mostrar solo las líneas de cabecera FASTA (>)
echo ">BRCA1_human
ATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC
>TP53_human
ATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA
>EGFR_human
ATGCGACCCTCCGGGACGGCCGGGGCAGCGCTCCTGGCGCTG" > genes.fasta

grep ">" genes.fasta

# Contar las secuencias (líneas de cabecera = número de secuencias)
grep -c ">" genes.fasta

La opción -c muestra únicamente el número de líneas que coinciden.

# Buscar un gen concreto
grep "BRCA1" genes.fasta

Comprensión del formato FASTQ

Los datos generados por los equipos de secuenciación reales se encuentran en formato FASTQ. A diferencia de FASTA, este formato incluye puntuaciones de calidad.

# Crear un archivo FASTQ (cada registro consta de cuatro líneas)
echo "@read_001
ATGCGATCGATCGATCGATCG
+
IIIIIIIIIIIIIIIIIIIII
@read_002
GCTAGCTAGCTAGCTAGCTAG
+
IIIIIIIIIIIIIIIIIIIII
@read_003
NNNNNATCGATCGATCGATCG
+
!!!!!IIIIIIIIIIIIIIIII" > reads.fastq

cat reads.fastq

Estructura de 4 líneas del archivo FASTQ:

1. Identificador de lectura que comienza con @
2. Secuencia de bases
3. + (línea separadora)
4. Puntuación de calidad (cuanto mayor, mejor; I = máximo, ! = mínimo)

awk: Procesamiento de datos

awk es una herramienta que procesa datos de texto columna por columna. Es similar a seleccionar columnas específicas en Excel.

# Crear datos separados por tabulaciones
echo "BRCA1	chr17	43044295	43170245
TP53	chr17	7661779	7687538
EGFR	chr7	55019017	55211628" > gene_locations.tsv

# Mostrar solo el nombre del gen (primera columna)
awk '{print $1}' gene_locations.tsv

# Calcular la longitud del gen (cuarta columna - tercera columna)
awk '{print $1, $4 - $3, "bp"}' gene_locations.tsv

Tubería (|): Conectar comandos

La tubería (|) conecta la salida de un comando con la entrada del siguiente, de forma similar a cómo se añadiría el resultado del reactivo A al reactivo B en un experimento.

# Extraer y contar solo las cabeceras del FASTA
echo ">BRCA1_human
ATGGATTTATCTGCTCTTCGCGTTGAAGAAGTACAAAATGTC
>TP53_human
ATGGAGGAGCCGCAGTCAGATCCTAGCGTGAGTTTGCTGTGA
>EGFR_human
ATGCGACCCTCCGGGACGGCCGGGGCAGCGCTCCTGGCGCTG" > genes.fasta

grep ">" genes.fasta | wc -l

Prueba por ti mismo (Ejemplo simplificado)

Completa el siguiente espacio en blanco para completar el comando que cuenta las secuencias en un archivo FASTQ.

En FASTQ, el identificador de secuencia comienza con @.

# ¿Cómo se cuentan las líneas que empiezan por @?
grep  "" test.fastq

Errores comunes y soluciones

P: Aparece command not found

En Colab, debes agregar el prefijo !: !grep ">" genes.fasta

P: El texto en coreano se muestra incorrectamente.

Es un problema de codificación de archivos. Verifica la codificación con file genes.fasta y, si no es UTF-8, conviértelo a UTF-8 usando iconv -f EUC-KR -t UTF-8 input.txt > output.txt.

P: Aparece Permission denied

Ocurre cuando el archivo no tiene permisos de ejecución. Asigna los permisos con chmod +x script.sh.

En el siguiente texto, aprenderás a analizar estos datos con mayor precisión utilizando Python.