La segunda columna de un archivo en formato SAM está reservado al denominado bitwise flag. Este es un número en formato decimal, que debe ser traducido a formato binario de 12 bits. Una vez traducido, cada uno de los bits, de derecha a izquierda, confirma (bit = 1) o rechaza (bit = 0) cada una de las declaraciones siguientes que se hacen sobre una lectura:
bit 1: La pareja ha tenido múltiples lecturas en la secuenciación
bit 2: Cada lectura del par ha sido correctamente alineada según el alineador
bit 3: Esta lectura no ha sido mapeada
bit 4: La pareja de esta lectura no ha sido mapeada
bit 5: La lectura mapea en la hebra reversa
bit 6: La pareja de la lectura mapea en la hebra reversa
bit 7: Esta lectura es la primera del par
bit 8: Esta lectura es la última del par
bit 9: Esta lectura tiene un alineamiento secundario
bit 10: Esta lectura no ha pasado los controles de calidad
bit 11: Esta lectura es un duplicado de PCR u óptico
bit 12: Esta lectura tiene un alineamiento suplementario
Este 'flag' puede ser muy útil. Fundamentalmente me ha servido hasta el momento para diferenciar la hebra de mapeo (bit #5) o para filtrar los duplicados de PCR (bit #11).
¿Como podemos acceder cómodamente a la información de este flag? Incluyo dos ejemplos de como seleccionar una lectura según su hebra de mapeo en Perl y Python que espero que sean útiles:
PERL
while (my $line = <SAM>) { # Lee el archivo SAM línea a línea
chomp $line;
if ($line =~ /^@/) {next;} # Descarta las líneas de cabecera
else {
my @chunks = split ("\t",$line);
my $bitwise = $chunks[1]; # Selecciona la columna del 'bitwise flag'
my $byte = sprintf("%012b", $bitwise); # Convierte el número en una cadena de 12 bits
my $strand = substr($byte, -5, 1); # Selecciona el bit que indica la hebra de la lectura
if ($strand) {print "La lectura mapea en la hebra reversa\n";}
else {print "La lectura mapea en la hebra directa\n";}
}
}
PYTHON
for line in sam_file:
if line[0] == '@': # Descarta las líneas de cabecera
continue
else:
chunks = line.split()
bitwise = int(chunks[1]) # Selecciona la columna del 'bitwise flag'
bits = bin(bitwise) # Transforma a binario, pero con un prefijo '0b'
bits = bits[2:] # Quita el prefijo '0b'
byte = "{0:0>12}".format(bits) # Añade ceros hasta obtener una cadena de 12 bits
strand = int(byte[-5]) # Selecciona el bit que indica la hebra de la lectura
if strand:
print "La lectura mapea en la hebra reversa"
else:
print "La lectura mapea en la hebra directa"
