Vous avez dit protocole ?
Qui dit biologie et bioinformatique dit protocole expérimental. C'est le cœur de la démarche scientifique, et un formalisme adapté est la clef pour assurer la reproductibilité des expériences, et ainsi garantir la validation des découvertes par la communauté. En paillasse, les solutions pour formaliser et conserver les protocoles sont plutôt naturellement pragmatiques et éprouvées par les années. Le papier a encore de beaux jours devant lui. Cependant, il est assez difficile d'adapter le principe du cahier de manip au monde de l'informatique, car il est inefficace tel quel.
À la place, beaucoup ont adopté les systèmes de gestion de workflow, qui permettent de retranscrire de manière assez naturelle la notion de reproductibilité au support informatique. Beaucoup d'outils spécialisés à la bioinformatique ont depuis vu le jour : Galaxy (dont j'avais déjà parlé ici), Taverna, Ergatis, etc. Ils permettent de créer, exécuter et partager des protocoles entiers, mais la plupart sont développés pour offrir une interface graphique à un utilisateur novice en programmation, et sont peu voire pas adaptés à l'utilisation en ligne de commande. De plus, l'utilisateur est souvent dépendant d'un ensemble d'outils pré-intégrés au système, ce qui est souvent handicapant.
Rendez-moi mon shell !
L'outil de référence pour créer des workflows en ligne de commande, vers lequel la plupart des gens se tournent encore aujourd'hui, est GNU Make. Il est bien connu, documenté, éprouvé, et permet de définir des workflows ramifiés. Un de ses points forts est le fait que si une partie du workflow seulement doit être modifiée, il n'exécutera que les parties nécessaires. Étant donnés les volumes des données manipulées en bioinformatique, c'est un avantage de taille.
Même s'il est généralement utilisé pour la compilation de logiciels, il est très adapté à la démarche scientifique. En effet, lorsque vous définissez un protocole expérimental, il y a de grandes chances pour qu'il soit à peu de choses près directement transposable en Makefile, les différentes étapes du protocole sont souvent assimilable à des règles pour GNU Make.
Mighty Python
Seulement voilà : GNU Make est un dinosaure poussiéreux sacrément compliqué à prendre en main, et long à maîtriser. Heureusement, il existe des initiatives dont le but est de proposer des alternatives plus modernes à cette honorable antiquité. C'est le cas, par exemple, de SCONS, Ant, Waf, Rake, Aap, Ninja et Snakemake. Même si je ne traiterai ici que de ce dernier, qui propose une solution similaire à GNU Make basé sur Python, il n'est pas exclu que d'autres alternatives puissent être aussi valables dans ce contexte.
Le principe est similaire : On définit des étapes/règles qui ont chacune des fichiers d'entrée, des fichiers de sortie, et des choses à exécuter pour produire les fichiers de sortie à partir des fichiers d'entrée. On les écrit dans un fichier nommé Snakefile (équivalent d'un Makefile chez GNU Make) et on lance l'exécution avec la commande snakemake. Il propose en plus de retranscrire le workflow en graphe acyclique orienté (DAG) pour mieux le visualiser, ce qui permet notamment de vérifier si le Snakefile correspond bien au protocole expérimental désiré.
Il n'est d'ailleurs tellement pas surprenant d'utiliser Snakemake dans ce contexte que le premier exemple fourni dans leur documentation est un exemple de bioinformatique. Ça ne s'invente pas. En voici la traduction :
Cufflinks est un outil permettant d'assembler des transcrits, estimer leur abondance, et effectuer des analyses d'expression différentielle sur des données RNA-Seq. Cet exemple montre comment créer un workflow typique pour Cufflinks avec snakemake. Il suppose que des données RNA-Seq alignées sont fournies sous format BAM pour quatre échantillons 101–104. Pour chaque échantillon, les transcrits sont assemblés avec '''cufflinks''' (rule assembly). Les assemblages sont fusionnés en un fichier GTF avec '''cuffmerge''' (rule merge_assemblies). Une comparaison avec la piste GTF de hg19 est ensuite effectuée (rule compare_assemblies). Enfin, les expressions différentielles sont calculées sur les transcrits trouvés (rule diffexp).
1 TRACK = 'hg19.gtf' # should be fetched from the cufflinks page since special fields are expected<br><br>REF = 'hg19.fa'<br><br>CLASS1 = '101 102'.split()<br>CLASS2 = '103 104'.split()<br>SAMPLES = CLASS1 + CLASS2<br><br>rule all :<br>input : 'diffexp/isoform_exp.diff', 'assembly/comparison'<br><br>rule assembly :<br>input : 'mapped/{sample}.bam'<br>output : 'assembly/{sample}/transcripts.gtf', dir='assembly/{sample}'<br>threads : 4<br>shell : 'cufflinks –num-threads {threads} ‑o {output.dir} –frag-bias-correct {REF} {input}'<br><br>rule compose_merge :<br>input : expand('assembly/{sample}/transcripts.gtf', sample=SAMPLES)<br>output : txt='assembly/assemblies.txt'<br>run :<br>with open(output.txt, 'w') as out :<br>out.write('\n'.join(input))<br><br>rule merge_assemblies :<br>input : 'assembly/assemblies.txt'<br>output : 'assembly/merged/merged.gtf', dir='assembly/merged'<br>shell : 'cuffmerge ‑o {output.dir} ‑s {REF} {input}'<br><br>rule compare_assemblies :<br>input : 'assembly/merged/merged.gtf'<br>output : 'assembly/comparison/all.stats', dir='assembly/comparison'<br>shell : 'cuffcompare ‑o {output.dir}all ‑s {REF} ‑r {TRACK} {input}'<br><br>rule diffexp :<br>input : expand('mapped/{sample}.bam', sample=SAMPLES), gtf='assembly/merged/merged.gtf'<br>output : 'diffexp/gene_exp.diff', 'diffexp/isoform_exp.diff'<br>threads : 8<br>run :<br>classes = [','.join(expand('mapped/{sample}.bam', sample=cls)) for cls in (CLASS1, CLASS2)]<br>shell('cuffdiff –num-threads {threads} {gtf} {classes[0]} {classes[1]}')On peut remarquer qu'il s'agit ni plus ni moins d'un script python : il est tout à fait possible d'utiliser tous les objets et fonctions disponibles dans le langage, et importer des modules en fonction des besoins. Pour chaque règle, il est possible de demander directement l'exécution de code Python directement intégré dans le Snakefile lorsqu'il est placé dans le bloc "run" (voir rule compose_merge). Bref, si vous êtes familier avec ce langage, vous ne serez pas dépaysé.
Retours d'expériences (pun intended)
Bien entendu, je ne vous en parlerais pas si je ne l'avais pas testé moi-même. Pour être tout à fait honnête, je pense qu'une solide expérience de développement en Python n'est pas du luxe. Je trouve Snakemake plus user-friendly que Make, mais il reste un système complexe qui demandera de la patience et quelques efforts à prendre en main. Cela dit, quand vous manipulez de grands volumes de données, c'est agréable de pouvoir facilement choisir quoi réexécuter en fonction des besoins.
Je l'utilise dans le cadre d'un projet de recherche où les protocoles ne sont pas pré-établis, et dans ce contexte l'avantage d'avoir un système facilement éditable et exécutable en ligne de commande se fait évident. J'apprécie beaucoup la souplesse du système, il faut avouer que c'est autre chose que les scripts shell cracras ou les Makefiles abscons. Snakemake a encore du chemin à faire, il y a pas mal de fonctionnalités qui manquent un peu ou semblent être encore en chantier, mais il est tout à fait utilisable en l'état, et je l'ai adopté.
L'option permettant de visualiser le pipeline sous forme de DAG est également très appréciable, car elle permet de visualiser l'état d'avancement de son exécution. Les étapes déjà effectuées s'affichent en pointillés et celles qui seront exécutées au prochain appel de la commande sont entourées d'un trait plein. Voici ce que ça donne pour un Snakemake de mon projet :
Je vous invite à visiter leur wiki si vous souhaitez en savoir plus. Snakemake est hébergé sur Bitbucket : https://bitbucket.org/johanneskoester/snakemake
______________________
Le premier à me dire ce qui ne va pas dans mon superbe protocole décrit dans la première illustration (à part les noms débiles) aura gagné le cocotier.
En effet, il est incompatible avec Snakemake. Pourquoi ?
______________________
Disclaimer : Tous les contenus utilisés pour les images finales de cet article sont soit dans le domaine public (1) soit sous licence Creative Commons Attribution-Share Alike 3.0 Unported (2,3,4), soit de ma production.
Merci à Sam pour m'avoir fait découvrir Snakemake et ainsi donné l'idée générale de cet article, et à mes relecteurs chéris pour leur patience, leurs correction et leurs idées : Guillaume, Mica, Nallias et ZaZoOo (sisi miss, au moins pour la patience!)
Laisser un commentaire