Ce document vous présente la démarche et la réflexion pour analyser un jeu de données de métabarcoding avec l’outil FROGS [1]. Les analyses ont été effectuées sur l’instance Galaxy de la plateforme Migale avec la version 4.1.0.

Avertissement

Les réponses aux questions de ce TP peuvent être visibles en appuyant sur le bouton Solution. À n’utiliser que si vous bloquez bien sûr !

1 Introduction à Galaxy

1.1 TP de découverte / rappels

Galaxy propose des tutoriels pour découvrir l’interface ou apprendre à utiliser certains outils. Ce tutoriel vous explique comment utiliser l’interface, importer des données et comment les organiser.

1.2 Instance Galaxy de Migale

Connectez-vous sur l’instance Galaxy de la plateforme Migale : galaxy.migale.inrae.fr. Votre identifiant et votre mot de passe vous ont été communiqués dans le mail d’activation de votre compte sur la plateforme.

Connectez-vous sur l’instance

Solution

Mise en garde

Vous devez maintenant être capable de :

vous connecter sur l’instance Galaxy de Migale
créer des historiques et les nommer
vous déplacer dans les historiques
uploader des fichiers présents en local ou sur le web

2 Contrôle qualité des données de séquençage

Créez l’historique 16S ALL FASTQ
Uploadez les 128 fichiers FASTQ issus de l’archive 16S.tar.gz disponible à l’adresse http://genome.jouy.inra.fr/formation/Metagenomics_03-2021/16S.tar.gz
Lancez FastQC [2] pour analyser la qualité des 64 échantillons (128 fichiers)

Avertissement

Il est possible de choisir plusieurs datasets en même temps pour leur appliquer la même action ou de créer une collection !

Solution

Lancez MultiQC [3] pour obtenir une rapport synthétique

Solution

Quel est le principal atout de cet outil ?

Solution

MultiQC permet d’aggréger une multitude de rapports en un seul. C’est très intéressant pour comparer les échantillons entre eux et constater très facilement s’ils sont semblables ou pas, en termes de nombre de reads, de qualité, de contenu…

Voilà un rapport sur de vraies données qui va vous permettre d’être confrontés à des informations réelles.

Est-ce que tous les reads sont de la même longueur ? Qu’est-ce que cela signifie ?

Solution

Non, les reads sont de longueur différente. Cela indique qu’un traitement a été effectué sur les données brutes. Il s’agit de comprendre pourquoi. Si vous les recevez ainsi, contactez la plateforme de séquençage pour obtenir des informations. Ici, les primers de séquençage ont été retirés.

La qualité des bases vous semble-t-elle correcte ?

Solution

Les données de séquençage Illumina ont généralement des profils équivalents, avec une qualité de base qui se dégrade vers la fin du read.

3 Gérer des données multiplexées

La plateforme de séquençage peut vous fournir les données de vos échantillons en un seul fichier FASTQ, associé à un fichier texte. Dans ce cas-là, c’est à vous de démultiplexer vos données pour arriver à un fichier FASTQ ↔︎ un échantillon

Créez un historique appelé MULTIPLEXED
Placez-y les fichiers suivants:
- http://genome.jouy.inra.fr/formation/Metagenomics_03-2021/Multiplexed/barcode.tabular
- http://genome.jouy.inra.fr/formation/Metagenomics_03-2021/Multiplexed/multiplex.fastq

Vous avez à votre disposition un fichier FASTQ contenant des reads appartenant à différents échantillons. La correspondance est indiquée dans le fichier barcode.tabular :

Sample	5’ seq	3’ seq
MgArd0001	ACAGCGT	TGTACGT
MgArd0009	ACAGTAG	TGTACGT
MgArd0017	ACGTCAG	TGTACGT
MgArd0029	ACTCAGT	TGTACGT
MgArd0038	ACTCGTC	TGTACGT
MgArd0046	AGCAGTC	TGTACGT
MgArd0054	AGCTATG	TGTACGT
MgArd0062	AGCTCGC	TGTACGT
MgArd0073	AGTATCT	TGTACGT
MgArd0081	AGTCTGC	TGTACGT

La 1ère colonne correspond au nom de l’échantillon, la seconde à la séquence de barcode sur l’extrémité 5’, la troisième à la séquence du barcode sur l’extrémité 3’.

Exemple :

AGCTATGACTGGGTGTAAGAGCTGTGATTGCTAACACTGTGGCCGGGCCAGGGCACCTGGATAAATCGGATTAGATACCCGGGTATGTACGT

À quel échantillon sera assignée cette séquence ?

Solution

MgArd0054

Lancez l’outil FROGS - demultiplex reads et choisir les paramètres adaptés

Solution

Barcode file : barcode.tabular Single or Paired-end reads : Single Select fastq dataset : multiplex.fastq Barcode mismatches : 0 Barcode on which end? : Both ends

Pourquoi n’y a-t-il pas 10 fichiers FASTQ en sortie ?

Solution

Pour organiser les données, les 10 fichiers FASTQ sont rangés dans une archive appelée demultiplexed.tar.gz. Si vous l’ouvrez, vous verrez que les 10 fichiers y sont.

Que deviennent les séquences dont les barcodes ne sont pas retrouvés ?

Solution

Elles sont stockées dans le dataset undemultiplexed.tar.gz qui contient les séquences non démultiplexées au format FASTQ.

Pourquoi n’y a-t-il aucune séquence ambigous ?

Solution

Car aucun mismatch n’a été toléré dans la séquence des barcodes. En autorisant des mismatches, on risque d’avoir des ambiguités entre plusieurs séquences, et le cas échéant l’outil ne sait pas décider entre les différents échantillons possibles. Ces séquences sont donc classées ambigous.

4 Préparer une archive pour FROGS

FROGS a besoin que tous les fichiers d’entrée lui soient fournis, mais contrairement aux collections, l’outil ne doît pas être lancé sur chaque fichier. Sans archive, il faut donc entrer tous les échantillons les uns après les autres. C’est fastidieux et peut causer des erreurs.

Une archive est une sorte de fichier qui contient plein de fichiers.

Cette archive doit avoir les propriétés suivantes :

Elle ne doit contenir que les fichiers FASTQ, pas de dossiers ni d’autres fichiers
Les fichiers FASTQ peuvent être compressés
Les fichiers FASTQ doivent être suffixés, si les données sont pairées, par *_R1.fastq.gz* et *_R2.fastq.gz*
Le nom de l’échantillon est construit à partir des chaines de caractères précédant ces suffixes.

En ligne de commande :

cd MyDir # les fichiers FASTQ sont dans ce répertoire
tar zcvf Archive.tar.gz *.fastq.gz

Ce document explique comment créer une archive sous Windows.

{{< pdf images/archive_creation_galaxy.pdf width=100% height="800px" >}}

Pour utiliser une archive dans FROGS preprocess il faut qu’elle soit uploadée au format tar (Datatype). Il suffit alors, dans le formulaire de choisir Archive dans la section Input type.

5 Analyse du jeu de données 16S

Il s’agit de résultats issus du travail de Chaillou et al., 2015 s’intéressant aux communautés bactériennes de produits de la mer et de produits carnés.

Les séquences que vous allez analyser sont des données simulées générées à partir des abondances estimées dans le papier.

Les 64 échantillons correspondent à 8 réplicats des 8 produits différents, générés sur la région 16S V1-V3.
Les primers sont 27F (5’-AGAGTTTGATCCTGGCTCAG-3’) 534R (5’-ATTACCGCGGCTGCTGG-3’)
Des reads de 2x300 bp ont été générés
Des erreurs ont été ajoutés dans les séquences suivant le modèle d’erreur Illumina classique
10% de séquences chimériques ont été ajoutées

Vous avez à dispostion un fichier de métadonnées (http://genome.jouy.inra.fr/formation/Metagenomics_03-2021/metadata_16S.tsv) :

et une archive contenant les 64x2=128 fichiers FASTQ compressés.

5.1 FROGS Preprocess

Créez l’historique 16S et uploadez l’archive que vous avez utilisée précédemment: 16S.tar.gz
Lancez FROGS preprocess. Utilisez les informations que nous vous avons fournies pour déterminer quels paramètres utiliser, ainsi que la documentation en bas du formulaire de soumission

Solution

Sequencer: Illumina
Input Type: Archive
Archive file: 16S.tar.gz
Reads already merged: No
Reads 1 size: 300
Reads 2 size: 300
mismatch rate: 0.1
merge software: vsearch
Would you like to keep unmerged reads?: No
Minimum amplicon size: 400
Maximum amplicon size: 580
Sequencing protocol: Illumina standard
5’ primer: AGAGTTTGATCCTGGCTCAG
3’ primer: CCAGCAGCCGCGGTAAT (!)

Explorez le rapport HTML généré
Quelle taille d’amplicons choisir ?

Solution

Quand vous ne savez pas, mettez des bornes basses et hautes puis allez regarder la distribution des séquences (bouton Display amplicon lengths) pour définir des bornes plus précises.

Lisez la documentation de l’outil (bas du formulaire) pour comprendre ce qui a été fait
Quelles informations sont présentes dans le fichier FASTA ?

Solution

Dans le fichier FASTA sont stockées les séquences des amplicons mais pas seulement

Exemple:

31_4222;size=30 AAAGCGCCGCTGAGGGATGGGCCTGCGTCCGATTAGCTTGTTGGTGGGGTAACGGCCCACCAAGGCGACGATCGGTAG…

31_4222 est le nom du read (dans la réalité il ressemble plus à M02765:8:000000000-AHM7P:1:1103:17002:15991) size=4 indique le nombre de séquences identiques trouvées dans tous les échantillons

Quelles informations sont présentes dans le fichier TSV ?

Solution

Ce fichier est la table d’abondance de chaque amplicon unique par échantillon

Exemple:

#id           BHT0.LOT01  BHT0.LOT3 
31_4222               10        20

L’amplicon est vu 10 fois dans l’échantillon BHT0.LOT01 et 20 fois dans l’autre échantillon.

Est-ce attendu que certaines paires de reads ne contiguent pas ?

Solution

Tout dépend de la taille attendue de l’amplicon qui a été séquencé. Généralement, les régions du 16S ciblées sont censées permettre un chevauchement des paires. Dans ce cas précis, les paires qui ne contiguent pas ne sont pas attendues et peuvent (doivent) être supprimées. Pour les ITS, il est essentiel de les conserver par contre car certains ITS peuvent dépasser 600 paires de bases. Si la taille attendue de votre amplicon est supérieure à 600 bp, vous ne deviez pas avoir de paires qui contiguent.

Qu’est-ce que la déréplication ?

Solution

La déréplication permet de garder un unique exemplaire de plusieurs séquences identiques. Le détail de l’abondance de cette séquence par échantillon est indiquée dans le fichier count.tsv

5.2 FROGS Clustering

Lancez l’outil FROGS clustering en spécidiant une valeur de d=1** et en ne réalisant pas l’étape de denoising**

Solution

Sequences file: FROGS preprocess FASTA file
Count file: FROGS preprocess TSV file
Aggregation distance: 1
Perform denoising clustering step: No

Quelles informations sont contenues dans les datasets générés ?

Solution

FASTA: séquences représentatives des clusters construits

Exemple:

Cluster_1 1:N:0:99 TACGTAGGTGGCAAGCGTTATCCGGAATTATTGGGCGTAAAGCG…

BIOM: fichier organisé de façon hiérarchique permettant de stocker entre autres les abondances des clusters par échantillon

ne pas chercher à le lire à l’oeil

TSV: une ligne correspond à un cluster. L’ensemble des identifiants des séquences qui le composent est indiqué (en pratique on ne s’en occupe presque jamais)

Qu’est-ce que l’étape de denoising ?

Solution

Le denoising consiste à effectuer un clustering à d=1 avant de faire un clustering avec un d plus important. Cela permet d’accélérer de façon significative le traitement sans trop changer le résultat.

Lancez l’outil FROGS clusters stat

Solution

Explorez le rapport HTML généré
Combien de séquences sont contenues dans le plus gros cluster ?

Solution

83,705

Combien de clusters sont composés d’une seule séquence ?

Solution

50,811

Quel pourcentage du total des clusters représentent-ils ?

Solution

98.39%

5.3 FROGS Remove chimera

Lancez l’outil FROGS remove chimera pour supprimer les éventuelles chimères que pourrait détecter vsearch [4].

Solution

Explorez le rapport HTML généré
Combien de clusters ont été supprimés ?

Solution

14,159 (27.4%)

Combien cela représente-t-il de séquences ?

Solution

14,609 (2.6%)

Qu’en concluez-vous ?

Solution

Les clusters chimériques sont majoritairement composés de peu de séquences.

Quelle est la plus grande abondance d’un cluster chimérique détecté ?

Solution

5.4 FROGS Cluster filters

Lancez l’outil FROGS filters pour ne garder que les clusters avec une abondance supérieure à 0.005% et présents dans au moins 4 échantillons. Utilisez également la détection des séquences phiX

Solution

Explorez le rapport HTML généré
Combien de clusters ont été supprimés et combien représentent-ils de séquences ?

Solution

98.7% des clusters sont supprimés ! Mais cela représente seulement 7% des séquences.

Quelle information est présente dans le fichier excluded.tsv ?

Solution

Chaque colonne représente un filtre et chaque ligne un cluster. Un cluster peut bien entendu se retrouver dans plusieurs colonnes.

Quel(s) est (sont) le(s) filtre(s) qui a (ont) permis de supprimer le plus de clusters ?

Solution

La grande majorité des clusters supprimés sont supprimés à la fois parce qu’ils sont vus dans moins de 4 échantillons et avec une abondance relative < 0.00005.

Ces échantillons étaient-ils contaminés par des séquences de phiX restantes ?

Solution

Non

Bonus / New functionality

**Utilisez l’information concernant les échantillons répliqués pour ne garder que les clusters présents dans au moins 50% d’un groupe. Le fichier est disponible ici : https://genoweb.toulouse.inra.fr/~formation/15_FROGS/current/chaillou_replicate_information.tsv

5.5 FROGS Taxonomic affiliation

Lancez l’outil FROGS taxonomic affiliation pour effectuer une affiliation Blast avec la base de données 16S SILVA Pintail100 138.1

Solution

Explorez le rapport HTML généré
Est-ce que tous les ASVs ont été affiliés à une séquence présente dans la base de référence utilisée ?

Solution

Non aucune !

Que feriez-vous des ASVs qui ne seraient pas affiliés ?

Solution

Il ne faut surtout pas les supprimer. Il s’agit peut-être de séquences non présentes dans la base de données utilisée. Dans ce cas, il faudrait utiliser d’autres bases de données (NCBI par exemple). Il est important de savoir ce qu’elles sont. Elles peuvent être une contamination (hôte, environnement, amplification d’autres espèces…). Cela donnera des informations importantes sur la qualité de votre expérimentation.

Expliquez le graphique Blast multi-affiliation summary

Solution

Ce graphique représente le nombre d’ASVs / séquences multi-affiliées à un rang donné. Cela signifie que plusieurs hits identiques de l’ASV sont ressortis avec une ambiguité au rang suivant. Par exemple, pour une multi-affiliation au rang Genus, il y a deux hits d’espèces différentes.

Le nombre de multi-affiliations au rang Species vous paraît-il surprenant ? Pourquoi ?

Solution

Non, le 16S ne permet généralement pas de discriminer les espèces au-delà du genre.

Lancez l’outil FROGS affiliation stats

Solution

Explorez le rapport HTML généré
Qu’est-ce qu’une courbe de raréfaction ? Quel est son intérêt principal ?

Solution

Une courbe de raréfaction permet d’évaluer si la profondeur de séquençage est suffisante. Cela permet de dire qu’on a suffisamment séquencé pour que la majorité des espèces de l’échantillon soient séquencées. Il faut dans l’idéal que la courbe atteigne un plateau.

À quoi vous ferait penser un ASV affilié avec un % de couverture inférieur à 80% ?

Solution

À une potentielle chimère !

Que conclure si les % d’identité de vos ASVs les plus abondants sont faibles ?

Solution

Vos séquences d’intérêt ne sont pas présentes dans la base de données utilisée.

5.6 FROGS Biom to TSV

Lancez l’outil FROGS biom to tsv en demandant l’extraction des mulit-affiliations et la présence des séquences dans le fichier de sortie

Solution

Abundance file: FROGS taxonomic affiliation BIOM Sequences file: FROGS filters FASTA Extract multi-alignments: Yes

Explorez les datasets générés
Combien de séquences possède le plus gros cluster ?

Solution

83,705 séquences

Combien de séquences possède le moins gros ASV ?

Solution

28 séquences

Pourquoi le Cluster_1 a plusieurs affiliations dans le fichier multi-affiliations sans pour autant être multi-affilié ?

Solution

Le Cluster_1 est affilié exclusivement à Brochotrix thermosphacta. Vous vous rappelez les copies de 16S dans les génomes bactériens ? Il n’y pas d’ambiguité sur l’espèce, donc pas de multi-affiliation ici.

Que pensez-vous des affiliations du Cluster_3 ?

Solution

Le Cluster_3 est quasi-exclusivement affilié à l’espèce Lactobacillus sakei. Seule une affiliation (Bacteria;Firmicutes;Bacilli;Lactobacillales;Lactobacillaceae;Lactobacillus;unknown species) est contraire. Cette affiliation pourrait donc être légitimement changée.

Pour vous faciliter la vie, nous avons développé une interface web : AffiliationExplorer qui vous permet de modifier certaines multi-affiliations si nécessaire. Pour cela, il vous faudra télécharger sur votre ordinateur le fichier BIOM et le fichier de multi-affiliations, ainsi que le FASTA si vous avez éventuellement besoin de la séquence pour faire vos modifications.

Solution

Demonstration pour modifier une affiliation existante

5.7 FROGS Tree

Lancez l’outil FROGS tree

Solution

Explorez le rapport HTML généré
Que faire si certains ASVs ne sont pas inclus dans l’arbre ?

Solution

Ceux-là seront inexploitables pour calculer des indices basés sur la distance phylogénétique.

6 Analyse du jeu de données ITS

Créez l’historique ITS et importez le fichier http://genome.jouy.inra.fr/formation/Metagenomics_03-2021/ITS1.tar.gz

Vous avez ici un jeu de données ITS de fromages (real) et de communautés synthétiques. Il contient une archive contenant les fichiers FASTQ et les métadonnées associées :

Ce sont des données issues du projet INRAE MEM METABARFOOD (PRJNA685292) dont l’objectif est d’évaluer la pertinence des marqueurs eucaryotes classiquement utilisés pour caractériser les communautés des écosystèmes alimentaires d’intérêt à INRAE.

Il s’agit ici d’échantillons de fromages réels (echantillon*), de MOCKs simples à complexes d’espèces d’intérêt
Les primers sont ITS1-F (5’-CTTGGTCATTTAGAGGAAGTAA-3’) ITS2 (5’-GCTGCGTTCTTCATCGATGC-3’)
Des reads de 2x250 bp ont été produits par Illumina Miseq

6.1 FROGS Preprocess

Lancez FROGS preprocess. Utilisez les informations que nous vous avons fournies pour déterminer quels paramètres utiliser, ainsi que la documentation en bas du formulaire de soumission
Quelle taille d’amplicons choisir ?

Solution

Pour les ITS, les tailles attendues vont de quelques dizaines de nucléotides à plusieurs kb ! Notre expérience nous a montré qu’entre 50 et 490 (ou 590 pour du 2x300), tous les ITS sont conservés.

Sequencer: Illumina Input Type: Archive Archive file: ITS1.tar.gz Reads already merged: No Reads 1 size: 250 Reads 2 size: 250 mismatch rate: 0.1 merge software: vsearch Would you like to keep unmerged reads?: Yes Minimum amplicon size: 50 Maximum amplicon size: 490 Sequencing protocol: Illumina standard 5’ primer: CTTGGTCATTTAGAGGAAGTAA 3’ primer: GCATCGATGAAGAACGCAGC (!)

Explorez le rapport HTML généré
Lancez ensuite FROGS clustering, FROGS clusters stats, FROGS remove chimera et FROGS filters avec la même réflexion que pour les données 16S.
Explorez les rapports HTML générés pour vérifier que tout s’est bien passé

6.2 FROGS ITSx

Lancez FROGS ITSx sans utiliser l’option permettant de tronquer les séquences de SSU et de 5.8S

Solution

Sequences file: FROGS filters FASTA Abundance file: FROGS filters BIOM ITS region: ITS1 Check only if sequence detected as ITS?: Yes

Explorez le rapport HTML généré

6.3 FROGS Taxonomic affiliation

Lancez FROGS Taxonomic affiliation en spécifiant la base de données ITS_UNITE_Fungi_8.0 et sans faire d’assignation RDP.

Solution

Reference database: ITS_UNITE_Fungi_8.0 Also perform RDP assignation: No Sequence file: FROGS ITSx FASTA Abundance file: FROGS ITSx BIOM

Explorez le rapport HTML généré.
Lancez ensuite FROGS Affiliations stats et FROGS Biom to TSV avec la même réflexion que pour les données 16S
Explorez le rapport HTML généré
Comment expliquez-vous les statistiques d’alignement représentées dans le graphique Alignment distribution du rapport HTML de FROGS Affiliations Stats ?

Solution

Le % de couverture est plus bas par rapport aux données 16S vues précédemment. Deux possibilités: - soit les régions flanquantes ne sont pas dans la base de données - soit la base de données est trop éloignée pour que les alignements soient bons

Comment jugez-vous l’affiliation de l’ASV majoritaire (Cluster_1) ?

Solution

67% de couverture est trop faible pour être “normal”. Dans ces cas-là, il faut utiliser une autre base de données pour rechercher la vraie taxonomie de cette séquence (si elle se trouve quelque part). Il faudra donc passer du temps à expertiser les affiliations “à la main”.

Est-ce que de longs ITS sont retrouvés ?

Solution

Non, ils auraient pu être reconnus grâce au tag “FROGS_combined” dans le nom de l’ASV.

7 Toy datasets

7.1 16S

Data from https://doi.org/10.1371/journal.pone.0204629
http://genome.jouy.inra.fr/~orue/data.tar.gz
16S V3-V4
V3F(5’-ACGGRAGGCWGCAGT-3’) and V4R (5’-TACCAGGGTATCTAATCCT-3’)

8 Références

1. Escudié F, Auer L, Bernard M, Mariadassou M, Cauquil L, Vidal K, et al. FROGS: Find, Rapidly, OTUs with Galaxy Solution. Bioinformatics. 2018;34:1287‑94. doi:10.1093/bioinformatics/btx791.

2. Andrews S. FastQC A Quality Control tool for High Throughput Sequence Data. http://www.bioinformatics.babraham.ac.uk/projects/fastqc/. 2010. http://www.bioinformatics.babraham.ac.uk/projects/fastqc/.

3. Ewels P, Magnusson M, Lundin S, Käller M. MultiQC: summarize analysis results for multiple tools and samples in a single report. Bioinformatics. 2016;32:3047‑8.

4. Rognes T, Flouri T, Nichols B, Quince C, Mahé F. VSEARCH: a versatile open source tool for metagenomics. PeerJ. 2016;4:e2584.

Réutilisation

CC BY-SA 4.0

A work by Migale Bioinformatics Facility
Université Paris-Saclay, INRAE, MaIAGE, 78350, Jouy-en-Josas, France
Université Paris-Saclay, INRAE, BioinfOmics, MIGALE bioinformatics facility, 78350, Jouy-en-Josas, France