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.

module load delly/0.8.3

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

delly | cat

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

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:

module load samtools/1.13

samtools faidx Zymoseptoria_tritici.fa

samtools index  mapping_illumina_chr10_500kb.bam

Et maintenant, utilisons delly pour le variant calling:

delly call -g Zymoseptoria_tritici.fa  -o SV_calling_illumina.bcf  mapping_illumina_chr10_500kb.bam

Si je veux regarder ce résultat:

head -n3 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

# Conversion en fichier vcf
module load bcftools/1.10.2

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:

grep "^#"  SV_calling_illumina.vcf | head

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

grep -v "^#"  SV_calling_illumina.vcf | wc -l

b) comptage du nombre de SVs

Combien de variants sont présents dans le fichier?

# Combien de variants ? 
grep -v -c "^#"  SV_calling_illumina.vcf

Combien de variants sont de bonne qualité?

# Combien de variants de bonne qualité ? 
grep -v "^#"  SV_calling_illumina.vcf| grep -v "LowQual" | wc -l

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

grep -v "^#"  SV_calling_illumina.vcf| grep -v "LowQual" | grep  "<DEL>"| wc -l

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

grep -v "^#"  SV_calling_illumina.vcf| grep -v "LowQual" | grep  "<DUP>"| head

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

grep -v "^#"  SV_calling_illumina.vcf | grep  "INV"

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

grep -v "^#"  SV_calling_illumina.vcf | grep  "BND"

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

grep -v "^#"  SV_calling_illumina.vcf |grep -v "LowQual" | grep  "<INS>"| wc -l

c) extraction des informations des délétions

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

grep -v "^#"  SV_calling_illumina.vcf |grep -v "LowQual" | grep  "<DEL>" > SV_call.good_del.vcf

Récupération du chromosome et de la position

#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

grep -v "^#" SV_calling_illumina.vcf |head -1

Récupération des autres informations

#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

# ici exemple sans le sed (et donc pas de tabulation à la place des ;)
grep -v "^#" SV_calling_illumina.vcf | grep -v "LowQual" | grep "<DEL>" | \
cut -f8 | cut -d ";" -f1,4,5,13 > delly_info_1.txt

Fusion des 2 fichiers et formattage:

#Fusion des deux fichiers
paste -d '\t' delly_del_start.txt delly_del_info.txt > delly_del.txt

#Formattage et ménage
 awk '{OFS="\t" ; print $1,$2,$4,$3,$5,$6}' delly_del.txt | sed "s/END=//g"  > delly_del.tsv

 awk '{OFS="\t" ; print $1,$2,$4,$3,$5,$6}' delly_del.txt | sed "s/END=//g"  > delly_del.tsv

Nettoyage:

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)

module load sniffles/1.0.11

sniffles --help

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

a) indexation du bam

samtools index mapping_minion_chr10_500kb.bam

b) Variant calling

sniffles -l 100 -m mapping_minion_chr10_500kb.bam -v SV_calling_minion.vcf

grep '^#' SV_calling_minion.vcf

c) comptage du nombre de SVs

Combien de SV sont positionées sur le chr10?

cat SV_calling_minion.vcf | grep ^chr_10 | wc -l

Combien sont des délétions?

cat SV_calling_minion.vcf | grep ^chr_10 | grep "DEL" | wc -l

Combien sont des duplications?

cat SV_calling_minion.vcf | grep ^chr_10 | grep "DUP" | wc -l

Combien d'entre elles sont des inversions?

cat SV_calling_minion.vcf | grep ^chr_10 | grep "INV" | wc -l

Combien d'entre elles sont des translocations?

cat SV_calling_minion.vcf | grep ^chr_10 | grep -w "TRA" | wc -l

Combien d'entre elles sont des insertions?

cat SV_calling_minion.vcf | grep ^chr_10 | grep "INS" | wc -l

d) extraction des positions des délétions

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

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

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

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

Coller tout ensemble

paste sniffles_del_start.txt sniffles_del_stop.txt sniffles_del_infos.txt > sniffles_del.tsv

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

cat delly_del.tsv

cat sniffles_del.tsv

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