Filtrage des reads selon leur longueur
Cette page est dédiée à comment filtrer et trier un fichier BAM pour ne sélectionner que les reads qui répondent à un critère spécifique de longueur.
Pour se faire nous aurons besoin d'outils de manipulation de données de séquençage comme samtools
ou awk
. Etant donné que nous avions déjà utilisé samtools
à plusieurs reprises je vous invite à consulter les pages précédentes pour le guide de l'installation. Passons cependant à l'installation de awk
.
Installation de awk :
Utilisation de samtools
et awk
samtools
et awk
Filtrer avec samtools
et awk
: Pour filtrer les reads par longueur, vous pouvez combiner samtools view
pour lire le fichier BAM et awk
pour filtrer les reads en fonction de leur longueur. Voici une commande exemple pour sélectionner les reads ayant une longueur comprise entre 1000 et 2000 bases :
Détails:
samtools view -h input.bam
: Cette commande lit le fichier BAMinput.bam
et inclut les en-têtes (-h
).awk
: Le scriptawk
filtre les lignes. Il passe les lignes d'en-tête (qui commencent par '@') et vérifie la longueur des séquences (la colonne$10
représente la séquence du read dans le format SAM).length($10) >= 1000 && length($10) <= 2000
: Cela spécifie que la longueur du read doit être d'au moins 1000 bases et au plus 2000 bases.samtools view -b -o filtered.bam
: Convertit la sortie en format BAM et écrit le résultat dansfiltered.bam
.
Optimisation avec des indices:
Pour améliorer l'efficacité, particulièrement avec de grands fichiers BAM, il est utile que le fichier BAM soit indexé (
input.bam.bai
). Vous pouvez indexer votre BAM avecsamtools index input.bam
avant de procéder au filtrage.
Validation:
Après le filtrage, il est recommandé de vérifier les statistiques du nouveau fichier BAM pour vous assurer que le filtrage a été effectué correctement, en utilisant
samtools flagstat filtered.bam
.
Pour filtrer les reads dans un fichier BAM de sorte que seuls les reads ayant une longueur supérieure à 11919 paires de bases (dans ccet exemple) soient conservés, vous pouvez utiliser la combinaison de samtools
et awk
pour effectuer ce filtrage spécifique. Voici comment procéder :
Commande pour filtrer les reads
Vous pouvez utiliser cette commande pour sélectionner uniquement les reads ayant une longueur supérieure à 11919 paires de bases :
Explication de la commande :
samtools view -h input.bam
: Affiche le contenu du fichier BAMinput.bam
en incluant les en-têtes (indiqué par l'option-h
).awk '{if($1 ~ /^@/ || length($10) > 11919) print}'
: Le scriptawk
passe à travers chaque ligne de la sortie desamtools
. Il laisse passer les en-têtes (les lignes qui commencent par@
). Pour les autres lignes, qui sont les reads, il vérifie si la longueur de la séquence du read ($10
représente la séquence du read dans le format SAM) est supérieure à 11919. Si c'est le cas, il imprime la ligne.samtools view -b -o filtered.bam
: Convertit la sortie filtrée en format BAM et l'écrit dansfiltered.bam
.
Vérification et indexation
Après avoir filtré les reads, il est recommandé de vérifier et d'indexer le fichier BAM résultant :
Vérifier le fichier filtré :
Cette commande affiche les premiers reads, s'il y en a, qui ne respectent pas le critère de longueur. Normalement, cette commande ne devrait rien retourner si le filtrage a été correctement effectué.
Indexer le fichier BAM filtré :
Ceci crée un index pour le fichier
filtered.bam
, ce qui est utile pour des accès rapides et efficaces lors des analyses subséquentes.
Ces étapes vous permettront de traiter efficacement votre fichier BAM en sélectionnant uniquement les reads de longueur supérieure à 11919 paires de bases, facilitant ainsi les analyses génomiques qui dépendent de longs reads.
Last updated