Quel qscore ?

Le Q-score est une métrique qui exprime la qualité d'une lecture de séquence individuelle. Le Q-score est généralement utilisé pour quantifier la confiance que l'on peut avoir dans chaque base appelée lors du séquençage.

Définition du Q-score

Q-score est une mesure logarithmique de la probabilité d'erreur de basecalling. Plus précisément, un Q-score est défini comme :

$$Q=−10log 10 ​ P$$

où 𝑃 est la probabilité qu'une base appelée soit incorrecte. Ainsi, un Q-score plus élevé indique une plus grande confiance dans la lecture de la base et vice versa. Par exemple, un Q-score de 20 indique une probabilité d'erreur de 1 sur 100, tandis qu'un Q-score de 30 indique une probabilité d'erreur de 1 sur 1000.

Utilisation du Q-score dans Dorado

Dans Dorado, comme dans d'autres outils de basecalling pour le séquençage par nanopores, le Q-score sert à évaluer et à filtrer les lectures en fonction de leur qualité. Dorado permet aux utilisateurs de spécifier un Q-score minimum (--min-qscore) lors du basecalling. Cela permet de filtrer et de ne conserver que les lectures qui atteignent ou dépassent un seuil de qualité spécifique, ce qui est utile pour éliminer les données de faible qualité qui pourraient compromettre les analyses en aval.

Importance du qscore

Filtrage des données

Le filtrage basé sur le Q-score aide à réduire le volume de données à analyser et améliore la précision des analyses en aval, telles que l'assemblage de génomes, l'identification de variants, et d'autres formes d'analyses bioinformatiques.

Comparaison de kits de séquençage

Le Q-score peut également être utilisé pour comparer les performances de différents kits de préparation d'échantillons ou de conditions expérimentales, offrant une mesure standardisée de la qualité des lectures.

Optimisation des protocoles

En surveillant les Q-scores, les chercheurs peuvent ajuster les protocoles de séquençage pour optimiser la qualité des données produites.

Analyse

Analyse du gène TTR à partir d'un amplicon de 1F à 4R (7096 paires de bases). C'est toutes les régions codantes et introns du gène TTR (Transthyretin). Ces résultats sont pris de 1 pod5 parmis les 48 (Mingon de nanopore (NATIF)) réalisés avec le modèle hac@5.0.0 et le kit SQK-NBD114-24. Pour obtenir ces différents résultats nous avons modifié le qscore au moment du basecalling. Toutes ces valeurs nous laissent à penser que le qscore réglé aux alentours des 17 permettent d'avoir des résultats fiables sans perte de profondeur. A voir si celà change sur tout le séquençage (les 48 pod5).

/home/grid/dorado-0.7.2-linux-x64/bin/dorado basecaller \
    -x "cuda:0" \
    --min-qscore 7 \
    --no-trim \
    --emit-fastq \
    /home/grid/dorado-0.7.2-linux-x64/bin/dna_r10.4.1_e8.2_400bps_hac@v5.0.0
    pod5/ | \
    /home/grid/dorado-0.7.2-linux-x64/bin/dorado demux \
    --kit-name SQK-NBD114-24 \
    --emit-fastq \
    --output-dir demultiplexed

DDN	N° dossier	Echantillon	FC	illumina	Nanopore	Profondeur (x) c.13	VAF(%) c.13 C>T	Profondeur (x) c.371	VAF(%)c.371 G>A		illumina	Nanopore	Profondeur (x) c.13	VAF(%) c.13 C>T	Profondeur (x) c.371	VAF(%) c.371 G>A
25/12/1945	2023-3726	RB16	FLG-114	RAS	RAS	1750	24	1258	11	Minion	RAS	RAS	16441	11	14670	11
05/11/1942	2024-0849	RB17	FLG-114	RAS	RAS	2683	24	1959	12	Minion	RAS	RAS	27768	12	26159	11
13/11/1962	2024-0851	RB18	FLG-114	POS/c.290C>A	POS/c.290C>A	1563	24	10681	12	Minion	POS/c.290C>A	POS/c.290C>A	11817	11	10202	11
05/01/1947	2024-0882	RB19	FLG-114	RAS	RAS	1313	26	1003	14	Minion	RAS	RAS	10599	11	9934	10
01/11/1948	2024-0869	RB20	FLG-114	RAS	RAS	2085	29	1749	12	Minion	RAS	RAS	23943	13	22449	11
18/08/1984	2024-0325	RB21	FLG-114	POS/c.290C>A	POS/c.290C>A	1556	23	1171	12	Minion	POS/c.290C>A	POS/c.290C>A	15205	10	13343	10
02/04/1970	220691295	RB22	FLG-114	POS/c.424G>A	POS/c.424G>A	1548	29	1374	13	Minion	POS/c.424G>A	POS/c.424G>A	15311	13	15079	11
16/03/1968	203301144	RB23	FLG-114	POS/c.424G>A	POS/c.424G>A	1595	25	1128	10	Minion	POS/c.424G>A	POS/c.424G>A	12567	11	11631	11

RB09 Gridion 1 seul bam (bam0)

qscore 0

Barcode	Nbr Bases	371 C > T %	Nbr Bases	13 C > T %	13 C > G %
1 POD
9	690	8	723	8	5
10	666	9	707	7	6
11	927	9	959	9	6

qscore 9 (par défaut)

Barcode	Nbr Bases	371 C > T %	Nbr Bases	13 C > T %	13 C > G %
1 POD
9	620	6	659	7	4
10	605	6	654	6	5
11	811	7	859	6	5
Tous POD
9	78 264	6	80 084	6	3
10	77 604	6	81 460	6	3
11	91 754	6	94 442	6	4

qscore 14

Barcode	Nbr Bases	371 C > T %	Nbr Bases	13 C > T %	13 C > G %
1 POD
9	562	5	589	5	3
10	539	6	585	4	4
11	711	6	770	4	5
Tous POD
9	69222	5	71956	5	3
10	68312	5	72931	5	3
11	80929	5	84629	5	3

qscore 15

Barcode	Nbr Bases	371 C > T %	Nbr Bases	13 C > T %	13 C > G %
1 POD
9	522	5	566	5	3
10	523	6	559	4	4
11	675	6	742	5	4

qscore 17

Barcode	Nbr Bases	371 C > T %	Nbr Bases	13 C > T %	13 C > G %
1 POD
9	436	5	465	5	3
10	444	5	455	4	4
11	567	5	619	4	3
Tous POD
9	55636	5	58835	4	2
10	54870	5	59291	4	2
11	65218	5	69171	4	3

RB10-Qscore17 bam merged

RB11-Qscore17 bam merged

qscore 19

Barcode	Nbr Bases	371 C > T %	Nbr Bases	13 C > T %	13 C > G %
1 POD
9	246	4	271	4	3
10	252	6	266	3	2
11	306	4	351	4	4

qscore 24

Il faut savoir que j'ai mis que 3 barcodes dans les tableaux du dessus mais en réalité j'en possède de 9 à 17 soit 9. Une fois le qscore passé à 24, j'ai perdu des barecodes et ils n'avaient presque pas de profondeur. Le qscore de 24 est donc trop élevé.

---

Pour effectuer des tests de qualité, nous pouvons automatiser tout ce processus via un script en .sh. Un point d'entrée de POD5 et on en ressort les BAM qui sont rangé en fonction du barcode et du qscore.

#!/bin/bash

# Définir le chemin de l'exécutable Dorado
DORADO_BIN="/home/grid/dorado-0.7.2-linux-x64/bin/dorado"
MODEL_PATH="/home/grid/dorado-0.7.2-linux-x64/bin/dna_r10.4.1_e8.2_400bps_hac@v5.0.0"

# Dossier contenant les données POD5
INPUT_DIR="pod5/"

# Scores de qualité à tester
QS_SCORES=(0 9 14 15 17 19 24)

# Boucle sur chaque score de qualité
for qscore in "${QS_SCORES[@]}"; do
    # Créer un répertoire pour les résultats de ce Q-score
    OUTPUT_DIR="demultiplexed_q${qscore}"
    mkdir -p "${OUTPUT_DIR}"

    # Basecalling et démultiplexage
    $DORADO_BIN basecaller -x "cuda:0" --min-qscore "$qscore" --no-trim --emit-fastq $MODEL_PATH $INPUT_DIR | \
    $DORADO_BIN demux --kit-name SQK-NBD114-24 --emit-fastq --output-dir "${OUTPUT_DIR}"

    # Dossier pour stocker les fichiers BAM
    BAM_DIR="bam_files_q${qscore}"
    mkdir -p "${BAM_DIR}"

    # Convertir chaque fichier FASTQ en BAM
    for fastq in "${OUTPUT_DIR}"/*.fastq; do
        bam_file="${BAM_DIR}/$(basename "${fastq}" .fastq).bam"
        minimap2 -ax map-ont -t 8 ref_genome.mmi "${fastq}" | samtools sort -o "${bam_file}"
    done
done

echo "Basecalling, démultiplexage, et conversion en BAM terminés pour tous les Q-scores."