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 :

sudo apt-get install gawk

Utilisation de 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 :

samtools view -h input.bam | awk '{if($1 ~ /^@/ || length($10) >= 1000 && length($10) <= 2000) print}' | samtools view -b -o filtered.bam
  • Détails:

    • samtools view -h input.bam : Cette commande lit le fichier BAM input.bam et inclut les en-têtes (-h).

    • awk : Le script awk 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 dans filtered.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 avec samtools 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 :

samtools view -h input.bam | awk '{if($1 ~ /^@/ || length($10) > 11919) print}' | samtools view -b -o filtered.bam

Explication de la commande :

  • samtools view -h input.bam : Affiche le contenu du fichier BAM input.bam en incluant les en-têtes (indiqué par l'option -h).

  • awk '{if($1 ~ /^@/ || length($10) > 11919) print}' : Le script awk passe à travers chaque ligne de la sortie de samtools. 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 dans filtered.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é :

    samtools view filtered.bam | awk '{if(length($10) <= 11919) print}' | head

    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é :

    samtools index filtered.bam

    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