Workflow - Structural Variant Calling¶

École de bioinformatique IFB-INSERM 2024

  • Gabriele Adam - INRAe

V - Variations structurales (SV) / TP¶

1. detection de SV à partir de données Illumina short-reads¶

En premier, nous allons utiliser Delly (papier ici doi: 10.1093/bioinformatics/bts378, code là: https://github.com/dellytools/delly). Pour commencer nous travailllons sur données short-reads.

In [ ]:
module load delly/0.8.3

Comment connaître la commande à utiliser pour lancer delly ?

In [ ]:
delly | cat

Et plus spécifiquement, comment lancer delly pour faire du variant calling ?

In [ ]:
delly call | cat

a) Indexation¶

Comme vu dans les cours précédents, nous allons commencerpar indexer notre alignement et notre référence.

Load samtools:

In [ ]:
module load samtools/1.13
In [ ]:
samtools faidx Zymoseptoria_tritici.fa
In [ ]:
samtools index mapping_illumina_chr10_500kb.bam

Et maintenant, utilisons delly pour le variant calling:

In [ ]:
delly call -g Zymoseptoria_tritici.fa  -o SV_calling_illumina.bcf  mapping_illumina_chr10_500kb.bam

Si je veux regarder ce résultat:

In [ ]:
head SV_calling_illumina.bcf

Ici le format n'est pas VCF mais BCF. Il faut donc le modifier.

SV_calling_illumina.bcf may be a binary file. See it anyway? n

In [ ]:
# Conversion en fichier vcf
module load bcftools/1.10.2
In [ ]:
bcftools view SV_calling_illumina.bcf > SV_calling_illumina.vcf

Et maintenant, nous pouvons visualiser nos variants.

Etape 1: Regarder seulement le header du fichier:

In [ ]:
grep "^#" SV_calling_illumina.vcf

Etape 2: regarder seulement les 3 premiers variants trouvés:

In [ ]:
grep -v "^#" SV_calling_illumina.vcf | head -n3

b) comptage du nombre de SVs¶

Combien de variants sont présents dans le fichier?

In [ ]:
# Combien de variants ? 
grep -v -c "^#" SV_calling_illumina.vcf

Combien de variants sont de bonne qualité?

In [ ]:
# Combien de variants de bonne qualité ? 
grep -v "^#" SV_calling_illumina.vcf | grep -v -c "LowQual"

Combien de variants de bonne qualité sont de type Deletion ?

In [ ]:
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep -c "<DEL>"

Combien de variants de bonne qualité sont de type Duplication ?

In [ ]:
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep -c "<DUP>"
In [ ]:
Combien de variants de bonne qualité sont de type Inversion?
In [ ]:
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep -c "<INV>"

Combien de variants de bonne qualité sont de type inter chromosomique translocation ?

In [ ]:
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep -c "<BND>"

Combien de variants de bonne qualité sont de type Insertion ?

In [ ]:
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep -c "<INS>"

c) extraction des informations des délétions¶

Extraire les délétions de bonne qualité:

In [ ]:
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep "<DEL>"

Récupération du chromosome et de la position

In [ ]:
#Récupération du start des variants
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep "<DEL>" | \
cut -f1,2 > delly_del_start.txt

Récupération des autres informations

In [ ]:
#Récupération des autres informations
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep "<DEL>" | \
cut -f8 | cut -d ";" -f1,4,5,13 | sed "s/;/\t/g" > delly_del_info.txt

Fusion des 2 fichiers et formattage:

In [ ]:
#Fusion des deux fichiers
paste -d '\t' delly_del_start.txt delly_del_info.txt > delly_del.txt
In [ ]:
#Formattage et ménage
 awk '{OFS="\t" ; print $1,$2,$4,$3,$5,$6}' delly_del.txt | sed "s/END=//g"  > delly_del.tsv

Nettoyage:

In [ ]:
rm delly_del_info.txt delly_del_start.txt delly_del.txt

2. detection de SV à partir de données MinION long-reads¶

Ici pour les données long reads, nous allons utiliser un autre caller: sniffles (papier ici https://www.nature.com/articles/s41592-018-0001-7, code là:https://github.com/fritzsedlazeck/Sniffles)

In [ ]:
module load sniffles/1.0.11
In [ ]:
sniffles --help

Pas besoin de génome de référence!

a) indexation du bam¶

In [ ]:
samtools index mapping_minion_chr10_500kb.bam

b) Variant calling¶

In [ ]:
sniffles -l 100 -m mapping_minion_chr10_500kb.bam -v SV_calling_minion.vcf
In [ ]:
grep '^#' SV_calling_minion.vcf

c) comptage du nombre de SVs¶

In [ ]:
Combien de SV sont positionées sur le chr10?
In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | wc -l

Combien sont des délétions?

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep "DEL" | wc -l

Combien sont des duplications?

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep "DUP" | wc -l

Combien d'entre elles sont des inversions?

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep "INV" | wc -l

Combien d'entre elles sont des translocations?

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep -w "TRA" | wc -l

Combien d'entre elles sont des insertions?

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep "INS" | wc -l

d) extraction des positions des délétions¶

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep "DEL" | head

Extraire le chromosome et la postion de début de la délétion

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep "DEL" | cut -f -1,2 > sniffles_del_start.txt

Extraire la position de fin de la délétion

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep "DEL" | cut -d ";" -f 4 | cut -d "=" -f 2 > sniffles_del_stop.txt

Extraire les informations complémentaires des délétions

In [ ]:
cat SV_calling_minion.vcf | grep ^chr_10 | grep "DEL" | cut -f 8 | cut -d ";" -f 1 > sniffles_del_infos.txt

Coller tout ensemble

In [ ]:
paste sniffles_del_start.txt sniffles_del_stop.txt sniffles_del_infos.txt > sniffles_del.tsv
In [ ]:
rm sniffles_del_start.txt sniffles_del_stop.txt sniffles_del_infos.txt

3. Comparaison des résultats de Delly et Sniffles¶

Ici nous allons observer les résultats de Delly et de Sniffles

In [ ]:
cat delly_del.tsv
In [ ]:
cat sniffles_del.tsv

Les délétions sont elles aux mêmes positions?