From 76e91af355e88a03155013694eead9c196f6a70c Mon Sep 17 00:00:00 2001 From: karinevilleneuve Date: Tue, 5 Mar 2024 10:04:01 -0500 Subject: [PATCH] add package ggplot2 --- docs/blast.html | 108 ++++++++++++++++++++--------------------- docs/search_index.json | 2 +- 2 files changed, 55 insertions(+), 55 deletions(-) diff --git a/docs/blast.html b/docs/blast.html index b5b3a2f..ceb6508 100644 --- a/docs/blast.html +++ b/docs/blast.html @@ -201,60 +201,60 @@

Chapitre 3 BLAST
install.packages("randomcoloR")
-install.packages("ggplot2")
-
-
-# Définir le nom de votre scientifique 
-ma_scientifique = "Abigail A. Salyers"
-
-# Importer le tableau de données
-tableau = read.table(file = "bio1410_2024.csv", sep=",", header=TRUE, check.names=FALSE) 
-# Définir le rang taxonomique d'intérêt 
-taxa_rank = "Genus"
-# Extraire vos échantillons 
-mes_échantillons = subset(tableau, Scientifique == ma_scientifique & Abundance > 0)
-
-# Afin de visualiser uniquementles 10 taxons les plus abondants on utilise ce qu'on appelle une *for loop* ou boucle en français 
-empt_df = list() # On génère une liste vide que l'on viendra populer avec les résultats de la foor loop
-i = 0 # Ensuite on définit un premier itérateur (de type numérique) qui compte le nombre de loop effectué
-# Début de la for loop
-for (région_échantillonnée in unique(mes_échantillons$Region)){ # Cette séquence est le début de la loop où l'itérateur de gauche (région) 
-    i = i + 1                                           # représente chaque item de l'itérateur de droite (mes_échantillons$Région)
-    sample = subset(mes_échantillons, Region == région_échantillonnée) # Générer une nouvelle table de données basée sur la variable définit par l'itérateur 
-    total_abundance = aggregate(sample$Abundance, by = list(taxa_rank = sample[[taxa_rank]]), FUN = sum) # Additionner ensemble les ASVs du même genres
-    top = head(total_abundance[order(total_abundance$x, decreasing = T),], n = 5) # filtrer en order croissant pour identifer les 10 plus abondants
-    others_df = sample[!sample[[taxa_rank]] %in% top$taxa_rank,] # identifier les autres genres ne faisant pas partie des 10 plus abonants 
-    others_list = others_df[[taxa_rank]] # extraire leur nom du tableau de données dans une liste  
-    sample[sample[[taxa_rank]] %in% others_list,][[taxa_rank]] = "Autres" # Renommer leur genre pour Autres
-    empt_df[[i]] = sample # Enregistrer ce nouveau tableau 
-}
-df = do.call("rbind",empt_df) # Combiner les deux tables de données en une seule 
-
-
-# Maintenant que nous avons notre tableau de données avec seulement les 10 taxons les plus abondants sous la langue et la main
-# nous pouvons commencer à travailler sur le graphique. On commence par générer une liste de couleurs associées à chacun des taxons
-
-n = nrow(data.frame(unique(df$Genus))) # On calcule le nobre de genre unique que nous avons 
-palette = randomcoloR::distinctColorPalette(n) # On génère une nouvelle palette 
-my_scale = ggplot2::scale_fill_manual(name = as.character(taxa_rank), values = palette, na.translate = FALSE, drop = TRUE, limits = force) # On assigne une couleur à chacun de nos genres uniques 
-
-# Nous pouvons finalement générer le graphique 
-graphique = ggplot2::ggplot(df, aes(x = Region, weight = Abundance, fill = .data[[taxa_rank]])) +
-  facet_grid(~ Region, scales = "free") +
-  geom_bar() +
-  labs(y = "Abondance relative (%)") +
-  scale_y_continuous(expand = c(0,0)) +
-  theme_classic() +
-  theme(text = element_text(size = 12), strip.text.x = element_blank(), 
-        plot.title = element_text(hjust = 0.5, margin = margin(0, 0, 50, 0))) + 
-  guides(fill = guide_legend(title = "Genre", title.position = "top")) + 
-  my_scale
-
-graphique # visualiser le graphique produit 
-
-# Enregistrer en format coma separeted value (csv) votre tableau d'abondance 
-write.csv(df, "df.csv", quote=FALSE)
+
install.packages("randomcoloR")
+install.packages("ggplot2")
+
+
+# Définir le nom de votre scientifique 
+ma_scientifique = "Abigail A. Salyers"
+
+# Importer le tableau de données
+tableau = read.table(file = "bio1410_2024.csv", sep=",", header=TRUE, check.names=FALSE) 
+# Définir le rang taxonomique d'intérêt 
+taxa_rank = "Genus"
+# Extraire vos échantillons 
+mes_échantillons = subset(tableau, Scientifique == ma_scientifique & Abundance > 0)
+
+# Afin de visualiser uniquementles 10 taxons les plus abondants on utilise ce qu'on appelle une *for loop* ou boucle en français 
+empt_df = list() # On génère une liste vide que l'on viendra populer avec les résultats de la foor loop
+i = 0 # Ensuite on définit un premier itérateur (de type numérique) qui compte le nombre de loop effectué
+# Début de la for loop
+for (région_échantillonnée in unique(mes_échantillons$Region)){ # Cette séquence est le début de la loop où l'itérateur de gauche (région) 
+    i = i + 1                                           # représente chaque item de l'itérateur de droite (mes_échantillons$Région)
+    sample = subset(mes_échantillons, Region == région_échantillonnée) # Générer une nouvelle table de données basée sur la variable définit par l'itérateur 
+    total_abundance = aggregate(sample$Abundance, by = list(taxa_rank = sample[[taxa_rank]]), FUN = sum) # Additionner ensemble les ASVs du même genres
+    top = head(total_abundance[order(total_abundance$x, decreasing = T),], n = 5) # filtrer en order croissant pour identifer les 10 plus abondants
+    others_df = sample[!sample[[taxa_rank]] %in% top$taxa_rank,] # identifier les autres genres ne faisant pas partie des 10 plus abonants 
+    others_list = others_df[[taxa_rank]] # extraire leur nom du tableau de données dans une liste  
+    sample[sample[[taxa_rank]] %in% others_list,][[taxa_rank]] = "Autres" # Renommer leur genre pour Autres
+    empt_df[[i]] = sample # Enregistrer ce nouveau tableau 
+}
+df = do.call("rbind",empt_df) # Combiner les deux tables de données en une seule 
+
+
+# Maintenant que nous avons notre tableau de données avec seulement les 10 taxons les plus abondants sous la langue et la main
+# nous pouvons commencer à travailler sur le graphique. On commence par générer une liste de couleurs associées à chacun des taxons
+
+n = nrow(data.frame(unique(df$Genus))) # On calcule le nobre de genre unique que nous avons 
+palette = randomcoloR::distinctColorPalette(n) # On génère une nouvelle palette 
+my_scale = ggplot2::scale_fill_manual(name = as.character(taxa_rank), values = palette, na.translate = FALSE, drop = TRUE, limits = force) # On assigne une couleur à chacun de nos genres uniques 
+
+# Nous pouvons finalement générer le graphique 
+graphique = ggplot2::ggplot(df, aes(x = Region, weight = Abundance, fill = .data[[taxa_rank]])) +
+  facet_grid(~ Region, scales = "free") +
+  geom_bar() +
+  labs(y = "Abondance relative (%)") +
+  scale_y_continuous(expand = c(0,0)) +
+  theme_classic() +
+  theme(text = element_text(size = 12), strip.text.x = element_blank(), 
+        plot.title = element_text(hjust = 0.5, margin = margin(0, 0, 50, 0))) + 
+  guides(fill = guide_legend(title = "Genre", title.position = "top")) + 
+  my_scale
+
+graphique # visualiser le graphique produit 
+
+# Enregistrer en format coma separeted value (csv) votre tableau d'abondance 
+write.csv(df, "df.csv", quote=FALSE)

Scientifiques

diff --git a/docs/search_index.json b/docs/search_index.json index 8ed393a..940a833 100644 --- a/docs/search_index.json +++ b/docs/search_index.json @@ -1 +1 @@ -[["blast.html", "Chapitre 3 BLAST", " Chapitre 3 BLAST Vous avez vue en classe une figure représentant l’abondance relative des principaux genres bactéries présents sous la langue ainsi que sur la paume de la main dominante de chacun des scientifiques de votre classe. Vous aimeriez maintenant générer une figure similaire à celle présentée en classe mais comprenant uniquement vos échantillons. Pour ce faire, vos démonstrateurs ont rédigé ce petit script R. Il vous suffit de suivre les étapes suivantes : Télécharger le tableau de données bio1410_2024.csv et l’enregistrer dans le dossier Documents de l’ordinateur que vous utilisez. Si vous utilisez votre ordinateur personnel plutôt que celui de l’école vos démonstrateurs pourrons vous aider à spécifier le dossier le travail. Créer un nouveau document de type R markdown. Insérer un nouveau bloc de code R et y copier/coller le script ci-dessous. À la première ligne du code, copier/coller le nom de votre scientifique à partir de la liste des scienfitiques disponible au bas de la page. Attention, assurez-vous de ne pas avoir d’espace entre les guillemets et le nom. Exécuter les commandes en appuyant sur le bouton vert en haut à droite du bloc de code ( ▶ ). install.packages("randomcoloR") install.packages("ggplot2") # Définir le nom de votre scientifique ma_scientifique = "Abigail A. Salyers" # Importer le tableau de données tableau = read.table(file = "bio1410_2024.csv", sep=",", header=TRUE, check.names=FALSE) # Définir le rang taxonomique d'intérêt taxa_rank = "Genus" # Extraire vos échantillons mes_échantillons = subset(tableau, Scientifique == ma_scientifique & Abundance > 0) # Afin de visualiser uniquementles 10 taxons les plus abondants on utilise ce qu'on appelle une *for loop* ou boucle en français empt_df = list() # On génère une liste vide que l'on viendra populer avec les résultats de la foor loop i = 0 # Ensuite on définit un premier itérateur (de type numérique) qui compte le nombre de loop effectué # Début de la for loop for (région_échantillonnée in unique(mes_échantillons$Region)){ # Cette séquence est le début de la loop où l'itérateur de gauche (région) i = i + 1 # représente chaque item de l'itérateur de droite (mes_échantillons$Région) sample = subset(mes_échantillons, Region == région_échantillonnée) # Générer une nouvelle table de données basée sur la variable définit par l'itérateur total_abundance = aggregate(sample$Abundance, by = list(taxa_rank = sample[[taxa_rank]]), FUN = sum) # Additionner ensemble les ASVs du même genres top = head(total_abundance[order(total_abundance$x, decreasing = T),], n = 5) # filtrer en order croissant pour identifer les 10 plus abondants others_df = sample[!sample[[taxa_rank]] %in% top$taxa_rank,] # identifier les autres genres ne faisant pas partie des 10 plus abonants others_list = others_df[[taxa_rank]] # extraire leur nom du tableau de données dans une liste sample[sample[[taxa_rank]] %in% others_list,][[taxa_rank]] = "Autres" # Renommer leur genre pour Autres empt_df[[i]] = sample # Enregistrer ce nouveau tableau } df = do.call("rbind",empt_df) # Combiner les deux tables de données en une seule # Maintenant que nous avons notre tableau de données avec seulement les 10 taxons les plus abondants sous la langue et la main # nous pouvons commencer à travailler sur le graphique. On commence par générer une liste de couleurs associées à chacun des taxons n = nrow(data.frame(unique(df$Genus))) # On calcule le nobre de genre unique que nous avons palette = randomcoloR::distinctColorPalette(n) # On génère une nouvelle palette my_scale = ggplot2::scale_fill_manual(name = as.character(taxa_rank), values = palette, na.translate = FALSE, drop = TRUE, limits = force) # On assigne une couleur à chacun de nos genres uniques # Nous pouvons finalement générer le graphique graphique = ggplot2::ggplot(df, aes(x = Region, weight = Abundance, fill = .data[[taxa_rank]])) + facet_grid(~ Region, scales = "free") + geom_bar() + labs(y = "Abondance relative (%)") + scale_y_continuous(expand = c(0,0)) + theme_classic() + theme(text = element_text(size = 12), strip.text.x = element_blank(), plot.title = element_text(hjust = 0.5, margin = margin(0, 0, 50, 0))) + guides(fill = guide_legend(title = "Genre", title.position = "top")) + my_scale graphique # visualiser le graphique produit # Enregistrer en format coma separeted value (csv) votre tableau d'abondance write.csv(df, "df.csv", quote=FALSE) Scientifiques A - J K - V Abigail A. Salyers Kalpana Chawla Alice Augusta Ball Khatijah Mohammad Yousoff Ana Roque de Duprey Leone Norwood Farrell Angelina Fanny Hesse Lise Meitner Augusta Ada Byron Lovelace Margaret Elaine Heafield Hamilton Chien-Shiung Wu Marie Maynard Daly Clark Creola Katherine Coleman Johnson Marjory Stephenson Dorothy Hodgkin Mary Golda Ross Dorothy Jean Johnson Vaughan Mary Jackson Elizabeth Bugie Gregory Nettie Stevens Elizabeth Garrett Anderson Rachel Carson Emmanuelle Marie Charpentier Rosalind Elsie Franklin Emmy Klieneberger-Nobel Roseli Ocampo-Friedmann Farah Alibay Ruth Ella Moore Floy Agnes Naranjo Stroud Lee Sally Kristen Ride Frances Joan Estelle Wagner Vera Florence Cooper Rubin Gladys Mae Brown West Jane Goodale Jessie Isabelle Price Johanna Westerdijk June Dalziel Hart Almeida "],["404.html", "Page not found", " Page not found The page you requested cannot be found (perhaps it was moved or renamed). You may want to try searching to find the page's new location, or use the table of contents to find the page you are looking for. "]] +[["index.html", "BIO1410 - 2024 - Analyses R Chapitre 1 Introduction", " BIO1410 - 2024 - Analyses R Karine Villeneuve Chapitre 1 Introduction Ce petit guide à été rédigé à l’intention des élèves du cours BIO1410 afin de familiariser les étudiants à l’analyse de données de séquençage d’amplicons avec R. Avant de débuter toute analyse, il est nécessaire de comprendre les principales fonctionnalités de R. La permière section Introduction à R Studio et R Markdown présente donc très brièvement comment utiliser R Mardown avec R Studio. Par la suite, la section BLAST comprend un petit script R permettant à l’élève d’identifier les 10 genres bactériens les plus abondants dans ses échantillons récoltés sous sa langue ainsi que sur la paume de sa mains. Ce script permet aussi de générer un tableau de données en format comma seperated value (csv) que l’étudiant devra utiliser pour obtenir les séquences d’acides nucléiques à rechercher sur la base de données du National Center for Biotechnology Information (NCBI) avec leur outil Basic Local Alignment Search Tool (BLAST). La section Traitement des fastqs présente ensuite comment les données brutes du séquençage d’amplicon du gène ARNs 16S ont été traitées dans R afin d’obtenir une matrice de l’abondance relative des variants de séquence d’amplicon (ASVs) par échantillon. Finalement, la section Statistiques présente l’ensemble des analyses statistiques réalisées avec R afin de comparer le microbiome des deux régions échantillonner. "],["introduction-à-r-studio-et-r-markdown.html", "Chapitre 2 Introduction à R Studio et R Markdown 2.1 Pourquoi utiliser R Markdown 2.2 Utliser R Markdowm", " Chapitre 2 Introduction à R Studio et R Markdown 2.1 Pourquoi utiliser R Markdown Plusieurs préfèrent rédiger leurs scripts dans des fichiers de type R Markdown plutôt que R Script pour diverses raions. Je recommande personnellement d’utiliser R Markdown pour les raisons suivantes : Ce type de fichier permet de facilement annoter son script entre les sections de codes qui sont comprises dans un bloc (“chunk” en anglais). Un script bien annoter permet non seulement au scripteur de s’y retrouver facilement mais aussi de partager son code à ses collègues ou même publiquement. En effet, il est de plus en plus commun de retrouver dans les articles scientifiques un lien vers les scripts générés et uilisés par les auteurs de l’article afin d’analyser leurs données. Mais aussi, les blocs de codes permettent d’exécuter seulement certaines sections de code à la fois, ce qui ultimement permet de mofidier puis exécuter seulement ces sections de code sans avoir à re-exécuter l’entièreté du script. 2.2 Utliser R Markdowm La première étape consiste donc à créer un nouveau document de type R Mardown. Pour cela, ouvrez RStudio, puis cliquez sur Fichier. Dans le menu déroulant sélectionner Nouveau Fichier puis R Markdown.... Dans la nouvelle fenêtre vous pouvez donner le titre que vous voulez à votre nouveau document. Par défault le nouveau document affiche une petite introduction ainsi que des exemples tel que sur l’image ci-dessous : Ces informations ne sont pas pertinentes et vous pouvez supprimer l’ensemble du texte sous l’entête (l’entête correspond à la section délimitée par les trois tirets ---). 2.2.1 Bloc de codes Dans R markdown, les lignes de codes à exécuter doivent être comprises dans un bloc de code. Le texte non compris dans un bloc n’est donc pas considéré comme du code, ce qui permet d’annoter minutieusement votre script entre les blocs afin de vous y retrouver facilement. Un bloc de code R doit toujours débuter avec les caractères suivants : ```{r} et se terminer avec les caractères suivant : ```. Un bloc de code ressemble donc à ceci : ```{r} ``` Un bloc de code peut être inséré avec l’une des façons suivantes : le raccourcit clavier : Ctrl + Alt + I tapper manuellement les caractères délimitants (```{r} ``` ) l’onglet Code puis Insert chunk le bouton vert avec le petit c et signe de plus en haut à droite. Un fois votre code rédigé dans le bloc, vous pouvez exécuter l’entièreté du code contenu dans ce block en appuyant sur le bouton vert en haut à droite du code ( ▶ ). Il est possible d’insérer des blocs de code de différents languages de programmation tels que Bash et Python, il suffit de remplacer le r entre les accolades par le nom du programme utilisé. Plusieurs autres options peuvent être appliqués sur les blocs, pour plus d’informations je vous recommande de consulter la documentation disponible sur internet. 2.2.2 Définir le répertoire de travail par défault Généralement, avant de commencer à analyser des données, il faut indiquer à R où se trouve les données sur notre ordinateur. Il existe plusieurs manières de procéder, dont notamment, spécifier au début du script un répertoire/dossier de travail par défault. Pour ce faire, une suffit d’utiliser la commande setwd() comme présenté sur l’exemple ci-dessous. Le texte entre les guillemets de l’exemple doit simplement être modifié pour le vrai chemin vers votre répertoire/dossier de travail. Une fois le repertoire de travail définit par défaut de la sorte il n’est pas nécessaire d’inclure le chemin vers les fichiers que l’on veut importer tant qu’il se trouve dans le répertoire spécifié. Cela s’applique aussi à la sauvegarde de tableaux, figures ou autres. ```{r} setwd("/chemin/vers/le/répertoire/") ``` 2.2.3 Importer des données Pour importer un fichier en format comma seperated value (csv) dans R on utilise la commande read.csv avec les arguments suivants : file : Spécifier le nom du fichier header : Indiquer si la première ligne de notre tableau correspond au nom des colonnes (valeurs possible TRUE/FALSE) row.names : Indiquer quelle colone de notre tableau correspond au nom des rangés sep : Indiquer quel caractère doit servir de séparateur pour les colonnes check.name : Ne pas systématiquement remplacer le trait d’union par un point (valeurs possible TRUE/FALSE) ```{r} df = read.csv(file = "nom_fichier.csv", header = TRUE, row.names = 1, sep = ",", check.names = FALSE) ``` Pour un tableau en format tab delimited on peut utiliser la fonction read.table et spécifier \\t comme séparateur. "],["blast.html", "Chapitre 3 BLAST", " Chapitre 3 BLAST Vous avez vue en classe une figure représentant l’abondance relative des principaux genres bactéries présents sous la langue ainsi que sur la paume de la main dominante de chacun des scientifiques de votre classe. Vous aimeriez maintenant générer une figure similaire à celle présentée en classe mais comprenant uniquement vos échantillons. Pour ce faire, vos démonstrateurs ont rédigé ce petit script R. Il vous suffit de suivre les étapes suivantes : Télécharger le tableau de données bio1410_2024.csv et l’enregistrer dans le dossier Documents de l’ordinateur que vous utilisez. Si vous utilisez votre ordinateur personnel plutôt que celui de l’école vos démonstrateurs pourrons vous aider à spécifier le dossier le travail. Créer un nouveau document de type R markdown. Insérer un nouveau bloc de code R et y copier/coller le script ci-dessous. À la première ligne du code, copier/coller le nom de votre scientifique à partir de la liste des scienfitiques disponible au bas de la page. Attention, assurez-vous de ne pas avoir d’espace entre les guillemets et le nom. Exécuter les commandes en appuyant sur le bouton vert en haut à droite du bloc de code ( ▶ ). install.packages("randomcoloR") install.packages("ggplot2") # Définir le nom de votre scientifique ma_scientifique = "Abigail A. Salyers" # Importer le tableau de données tableau = read.table(file = "bio1410_2024.csv", sep=",", header=TRUE, check.names=FALSE) # Définir le rang taxonomique d'intérêt taxa_rank = "Genus" # Extraire vos échantillons mes_échantillons = subset(tableau, Scientifique == ma_scientifique & Abundance > 0) # Afin de visualiser uniquementles 10 taxons les plus abondants on utilise ce qu'on appelle une *for loop* ou boucle en français empt_df = list() # On génère une liste vide que l'on viendra populer avec les résultats de la foor loop i = 0 # Ensuite on définit un premier itérateur (de type numérique) qui compte le nombre de loop effectué # Début de la for loop for (région_échantillonnée in unique(mes_échantillons$Region)){ # Cette séquence est le début de la loop où l'itérateur de gauche (région) i = i + 1 # représente chaque item de l'itérateur de droite (mes_échantillons$Région) sample = subset(mes_échantillons, Region == région_échantillonnée) # Générer une nouvelle table de données basée sur la variable définit par l'itérateur total_abundance = aggregate(sample$Abundance, by = list(taxa_rank = sample[[taxa_rank]]), FUN = sum) # Additionner ensemble les ASVs du même genres top = head(total_abundance[order(total_abundance$x, decreasing = T),], n = 5) # filtrer en order croissant pour identifer les 10 plus abondants others_df = sample[!sample[[taxa_rank]] %in% top$taxa_rank,] # identifier les autres genres ne faisant pas partie des 10 plus abonants others_list = others_df[[taxa_rank]] # extraire leur nom du tableau de données dans une liste sample[sample[[taxa_rank]] %in% others_list,][[taxa_rank]] = "Autres" # Renommer leur genre pour Autres empt_df[[i]] = sample # Enregistrer ce nouveau tableau } df = do.call("rbind",empt_df) # Combiner les deux tables de données en une seule # Maintenant que nous avons notre tableau de données avec seulement les 10 taxons les plus abondants sous la langue et la main # nous pouvons commencer à travailler sur le graphique. On commence par générer une liste de couleurs associées à chacun des taxons n = nrow(data.frame(unique(df$Genus))) # On calcule le nobre de genre unique que nous avons palette = randomcoloR::distinctColorPalette(n) # On génère une nouvelle palette my_scale = ggplot2::scale_fill_manual(name = as.character(taxa_rank), values = palette, na.translate = FALSE, drop = TRUE, limits = force) # On assigne une couleur à chacun de nos genres uniques # Nous pouvons finalement générer le graphique graphique = ggplot2::ggplot(df, aes(x = Region, weight = Abundance, fill = .data[[taxa_rank]])) + facet_grid(~ Region, scales = "free") + geom_bar() + labs(y = "Abondance relative (%)") + scale_y_continuous(expand = c(0,0)) + theme_classic() + theme(text = element_text(size = 12), strip.text.x = element_blank(), plot.title = element_text(hjust = 0.5, margin = margin(0, 0, 50, 0))) + guides(fill = guide_legend(title = "Genre", title.position = "top")) + my_scale graphique # visualiser le graphique produit # Enregistrer en format coma separeted value (csv) votre tableau d'abondance write.csv(df, "df.csv", quote=FALSE) Scientifiques A - J K - V Abigail A. Salyers Kalpana Chawla Alice Augusta Ball Khatijah Mohammad Yousoff Ana Roque de Duprey Leone Norwood Farrell Angelina Fanny Hesse Lise Meitner Augusta Ada Byron Lovelace Margaret Elaine Heafield Hamilton Chien-Shiung Wu Marie Maynard Daly Clark Creola Katherine Coleman Johnson Marjory Stephenson Dorothy Hodgkin Mary Golda Ross Dorothy Jean Johnson Vaughan Mary Jackson Elizabeth Bugie Gregory Nettie Stevens Elizabeth Garrett Anderson Rachel Carson Emmanuelle Marie Charpentier Rosalind Elsie Franklin Emmy Klieneberger-Nobel Roseli Ocampo-Friedmann Farah Alibay Ruth Ella Moore Floy Agnes Naranjo Stroud Lee Sally Kristen Ride Frances Joan Estelle Wagner Vera Florence Cooper Rubin Gladys Mae Brown West Jane Goodale Jessie Isabelle Price Johanna Westerdijk June Dalziel Hart Almeida "],["traitement-des-fastqs.html", "Chapitre 4 Traitement des fastqs 4.1 Pré-analyse 4.2 Inspecter la qualité, filtrer et rogner 4.3 Retirer les erreurs de nucléotides 4.4 Gnérer la matrice d’abondance d’ASVs par échantillon 4.5 Classification taxonomique 4.6 Obtenir les séquences d’acides nucléiques (ADN) 4.7 Sauvegarder les résultats en format comma seperated value (csv) 4.8 Raréfier 4.9 Sauvegarder les résultats en format comma seperated value (csv)", " Chapitre 4 Traitement des fastqs DADA2 est un pipeline bio-informatique qui comprend une série d’étapes permettant de filtrer les séquences brutes obtenues par le séquençage d’amplicons. Ultimement, le pipeline permet d’obtenir une matrice d’abondance d’ASVs en fonction des échantillons ainsi que la classification taxonomique des séquences associées aux différents ASVs. DADA2 a deux particularités qui le distingue des autres pipeline courament utilisés. D’une part, il procède la à modélisation de l’erreur dûe au séquençage ce qui est censé permettre des distinguer les séquences mutantes, des séquences érronées. D’autre part, contrairement à d’autres pipelines comme QIIME ou Mothur, DADA2 ne regroupe pas les séquences similaires à 97% en Unités Taxonomiques Opérationelles (OTUs). Ses ASVs ne subissent pas de regroupement si les séqences ne sont pas identiques à 100%. 4.1 Pré-analyse Importer les librairies que nous utilisons pour traiter les données library(tidyr) library(dplyr) library(dada2) library(phyloseq) library(DECIPHER) Identifier les échantillons présents dans le dossier de travail et extraire le nom de ces deniers sans les artefacts ajoutés lors du séquençage. # Créer une variable qui définit le lien vers le dossier contenant tous les fichiers en format fastq path = "~/BIO1410/2023/raw_fastqs" # Identifier tous les fichiers fastq du brin Forward fnFs = sort(list.files(path, pattern="_R1_001.fastq", full.names = TRUE)) # Identifier tous les fichiers fastq du brin Reverse fnRs = sort(list.files(path, pattern="_R2_001.fastq", full.names = TRUE)) # Extraire le nom de l'échantillon à partir du nom du fichier en assumant # que le format du nom de l'échantillon est le suivant : NOM-ÉCHANTILLON_XXX.fastq.gz sample.names = sapply(strsplit(basename(fnFs), "_"), `[`, 1) # Appliquer ce nom sans les artefacts du séquençage à nos variables names(fnFs) = sample.names names(fnRs) = sample.names 4.2 Inspecter la qualité, filtrer et rogner Inspecter la qualités des séquences Afin d’inspecter la qualité des séquences nous utilisons la commande plotQualityProfile. Cette commande représente graphiquement le Q score (axe des Y) associé à chaque nucléotide (axe des X). La médianne est en vert, et les quartiles en orange pointillé. Le Q score est une mesure qui nous renseigne sur la précision du séquençage (voir tableau ci-dessous). Q score Précision 10 90 % 20 99 % 30 99.9 % 40 99.9 % plotQualityProfile(fnFs, aggregate = TRUE) plotQualityProfile(fnRs, aggregate = TRUE) L’analyse de ces graphiques nous permet donc de choisir les paramètres de filtrage et de rognage de la prochaine étape. La qualité des séquences du brin Forward est généralement de meilleur qualité que celle des séquences du brin Reverse. L’intégralité du brin Forward (250 nucléotides) est d’assez bonnes qualité pour le conserver entier alors que pour le brin Reverse nous supprimerons les nucléotides après la position 240. En observant les graphiques générés, vous pouvesz constater que le début des séquences est aussi de moins bonne qualité et comprend les amorces que nous devons retirer. Filtrer et rogner les séquences Tout d’abord nous créons un dossier filtered ainsi que les objets filtFs filtRs pour stoquer les séquences filtrées. filtFs = file.path(path, "filtered", paste0(sample.names, "_F_filt.fastq.gz")) filtRs = file.path(path, "filtered", paste0(sample.names, "_R_filt.fastq.gz")) names(filtFs) = sample.names names(filtRs) = sample.names Nous utilisons la fonction filterAndTrim afin de filtrer et rogner les séquences. truncQ : définit un l’indice Q score minimale. A la première instance d’un score de qualité inférieur ou égal à truncQ, la séquence est tronquée. truncLen : définit à quelle longueur les séquences vont être tronquées. Les séquences plus courtes que la longueur sont éliminées. trimLeft : définit la longueur que l’on coupe du côté 5’ des séquences. Cela permet d’enlever les amorces si elles n’ont pas été préalablement retiré. 799F: AACMGGATTAGATACCCKG = 19 nucléotides / 1115R: AGGGTTGCGCTCGTTG = 19 nucléotides maxEE : définit le nombre maximum d’erreurs attendues autorisées dans une lecture. Ce filtre se base sur l’indice Q score. Plus on augmente le chiffre, moins on est strict. # Filtrer les séquences out = filterAndTrim(fnFs, filtFs, fnRs, filtRs, trimLeft = c(19,16), truncLen=c(250,240), maxN=0, maxEE=2, truncQ=2,rm.phix=TRUE, compress=TRUE, multithread=TRUE) 4.3 Retirer les erreurs de nucléotides Cet étape consiste à estimer le taux d’erreur de séquençage. Son but est de permettre de différencier les séquences mutantes et les séquences érronées.Le modèle d’erreur est calculé en alternant l’estimation des taux d’erreur et l’inférence de la composition de l’échantillon jusqu’à ce qu’ils convergent vers une solution cohérente. errF = learnErrors(filtFs, multithread=TRUE, randomize=TRUE) errR = learnErrors(filtRs, multithread=TRUE, randomize=TRUE) Nous pouvons maintenant appliquer l’algorithme entraîné sur nos séquences afin de retirer les nulcéotides erronés dadaFs = dada(filtFs, err = errF, pool = "pseudo", multithread=TRUE) dadaRs = dada(filtRs, err = errR, pool = "pseudo", multithread=TRUE) 4.4 Gnérer la matrice d’abondance d’ASVs par échantillon Tout l’intérêt du séquençage en paire réside dans le but de fusionner les deux brins afin d’accroitre notre confiance en leur fiabilité. La fusion permet également d’obtenir des amplicons plus long. mergers = mergePairs(dadaFs, filtFs, dadaRs, filtRs) Nous avons enfin notre matrice d’abondance d’ASVs que nous allons stoquer dans l’objet seqtab. # Construire la table d'ASV avec les séquences combinées seqtab = makeSequenceTable(mergers) L’étape de retrait des chimère vise à éliminer toutes les séquences non biologiques qui sont produites lorsqu’une séquence n’est pas complètement amplifié lors de la PCR et qu’elle s’hybride avec l’amorce lors du prochain cycle. seqtab.nochim = removeBimeraDenovo(seqtab, method = "consensus", multithread = TRUE) # Comparer la dimension de la tables d'ASV avec et sans chimières dim(seqtab.nochim) sum(seqtab.nochim)/sum(seqtab) 4.5 Classification taxonomique Grâce à l’implémentation de la méthode de classification naïve bayésienne, la fonction assignTaxonomy prend en entrée l’ensemble de séquences à classer ainsi qu’un ensemble de séquences de référence avec une taxonomie connue. Nous utilisons la base de données d’amplicons de Silva. taxa = assignTaxonomy(seqtab.nochim, "~/16S_Project/silva_nr99_v138.1_train_set.fa.gz", multithread=TRUE, tryRC=TRUE) Comme nous avons plusieurs ASVs non classifié à différents rangs taxonomique, nous allons venir copier le rang taxonomique supérieur suivant du mot “Unclassified” aux ASVs non-classifiés. taxid = data.frame(t(taxa)) taxid[] = lapply(taxid, as.character) taxa2 = tidyr::fill(taxid, colnames(taxid),.direction = "down") taxa2 = sapply(taxa2, function(x){paste0("Unclassified_", x)}) taxid[is.na(taxid)] = taxa2[is.na(taxid)] taxid = t(taxid) # Retirer les ASVs n'appartenant pas au domaine des Bactéries taxid = subset(as.data.frame(taxid), Kingdom =="Bacteria") seqtab.nochim = seqtab.nochim[,colnames(seqtab.nochim) %in% rownames(taxid)] 4.6 Obtenir les séquences d’acides nucléiques (ADN) Nous voulons conserver la séquences d’ADN correspondant à chacun des ASVs afin de pouvoir ultimement utiliser cette séquences avec la base de données BLAST de NCBI. # Importer le tableau de données des métadonnées meta = read.table("/home/kvilleneuve/16S_Project/BIO1410/metadata.csv",sep=",", row.names=1, header=TRUE) # Combiner la matrice d'aboncance, le tableau de métadonnées et la table de taxonomique dans un object de type phyloseq ps = phyloseq(otu_table(t(seqtab.nochim), taxa_are_rows=TRUE), tax_table(as.matrix(taxid)), sample_data(meta)) # On dna = Biostrings::DNAStringSet(taxa_names(ps)) names(dna) = taxa_names(ps) ps1 = merge_phyloseq(ps, dna) taxa_names(ps1) = paste0("ASV", seq(ntaxa(ps1))) 4.7 Sauvegarder les résultats en format comma seperated value (csv) # La matrice d'anbondance d'ASV par échantillon write.csv(seqtab.nochim, file = "~/BIO1410/raw_asv.csv") # Le tableau de classification taxonomique des ASVs write.csv(taxid, file = "~/BIO1410/raw_sequences.csv") # Les séquences d'ADN associées à chacun des ASVs ps1 %>% refseq() %>% Biostrings::writeXStringSet("~/BIO1410/raw_refseq.fna", append=FALSE, compress=FALSE, compression_level=NA, format="fasta") 4.8 Raréfier Retirer les taxons peu abondant Une des raisons pour retirer les taxons peu abondant est de limiter les temps passé à analyser des taxons rarement vus. Cela s’avère généralement une bonne pratique afin de filtrer le “bruit de fond” inutile (taxons qui ne sont en fait que des artefacts du processus de collecte de données). Nous élminons ici les taxons dont l’abondance relative est inférieur à 0.005 % tel que recommandé par Bokulich et al., 2013 # Define threshold for low abundant taxa minTotRelAbun = 5e-5 # Get total sum for each taxa x = taxa_sums(ps) # Identify taxa with a total sum greater than the defined threshold keepTaxa = (x / sum(x)) > minTotRelAbun # Filter out from the phyloseq object any taxa not identified in the keepTaxa object prunedSet = prune_taxa(keepTaxa, ps) # View how many taxa were removed by sample loss_taxa=data.frame(sample_sums(prunedSet), sample_sums(ps), (sample_sums(prunedSet)-sample_sums(ps))) Nous raréfions ensuite au premier nombre de séquence au-dessus de 1000. Dans note cas cela correspond à 1118. # Get dataframe of sequences per sample sample_size = as.data.frame(sample_sums(prunedSet)) # Filter for the lowest number above 1000 rare_value = sample_size[which.max((sample_size[,1] >= 1000)/sample_size[,1]),] # Rarefy to value identified as rare_value ps_rare = rarefy_even_depth(prunedSet, rare_value, rngseed = 112, replace = FALSE, trimOTUs = TRUE, verbose = TRUE) # Confirm rarefaction as a sanity check sample_sums(ps_rare) 4.9 Sauvegarder les résultats en format comma seperated value (csv) write.csv(as.data.frame(as(tax_table(ps_rare), "matrix")), file = "~/rarefied_taxa.csv") write.csv(as.data.frame(as(otu_table(ps_rare), "matrix")),file = "~/rarefied_asv.csv") "],["statistiques.html", "Chapitre 5 Statistiques 5.1 Thème pour les graphiques 5.2 Indice de Shannon 5.3 Ordination (PCoA) 5.4 LEfSe 5.5 Abondances relatives", " Chapitre 5 Statistiques Importer les libraires que nous devons utiliser library(vegan) library(tidyverse) library(devtools) library(forcats) # To reorder our factors (function "fct_relevel") library(dplyr) # Dataframe manipulation (notably function "pull") library(tidyr) # Dataframe manipulation (function "separate") library(phyloseq) # Very pratical library for the analysis of amplicon data library(randomcoloR) # Generate sets of random colors library(ggplot2) # Generate plots library(stringr) # Makes working with strings as easy as possible (function "str_replace") library(ggtext) # Allows the use of markdown text (used to make names italic) library(ggpubr) 5.1 Thème pour les graphiques Je génère mon propre thème pour les différents graphiques à générer dans lequel je spécifie par exemple la taille du texte, la taille des lignes des axes, etc… text_size = 12 custom_theme = function(){ theme_classic() %+replace% theme( #text elements plot.title = element_text(size = text_size), #set font size plot.subtitle = element_text(size = text_size), #font size plot.caption = element_text(size = text_size), #right align axis.title = element_text(size = text_size), #font size axis.text = element_text(size = text_size), #font size axis.text.x = element_text(size = text_size), axis.text.y = element_text(size = text_size), legend.text = element_text(size = text_size), legend.title = element_text(size = text_size, hjust=0), strip.text = element_text(size = text_size), strip.background = element_blank(), panel.border = element_rect(colour = "black", fill = NA, linewidth = 1), panel.grid.minor = element_blank(), panel.grid.major = element_line(linewidth = 0.25, linetype = 'solid', colour = "grey") ) } Importer les fichiers produit par le pipeline DADA2 (incluant la raréfaction) asv = read.table(file = "~/rarefied_asv.csv", sep=",", row.names=1, header=TRUE, check.names=FALSE) taxa = read.table(file = "~/rarefied_taxa.csv", sep=",", row.names=1, header=TRUE) meta = read.table(file = "~/metadata.csv", sep=",", row.names=1, header=TRUE) # Insérer le symbole "\\n" (signifie saut de line) afin que le texte s'affiche sur deux lignes dans la légende de notre figure meta$Region = gsub("Paume de la main dominante","Paume de la \\nmain dominante", meta$Region) # Combiner les tableaux dans un objet de type phyloseq ps = phyloseq(otu_table(asv, taxa_are_rows = TRUE), tax_table(as.matrix(taxa)), sample_data(meta)) 5.2 Indice de Shannon Calculer l’indice de diversité par échantillon shannon = cbind(estimate_richness(ps, measures = 'shannon'), sample_data(ps)) Maintenant que nous avons les valeurs de diversité (indice de Shannon) nous aimerions comparer ces valeurs en fonction de la région échantilloné (langue et mains). Pour ce faire, nous pouvons utiliser un test de t de student, mais nous devons avant nous assurer que nos données respectent les postulats de ce test (distribution normale et variance similaire) # Est ce que les données ont une distibution gausienne (normale) ? valeur_shapiro=list() # créer une liste vide # Réaliser le test de shapiro pour chacune des régions for (region in unique(shannon$Region)){ sample = subset(shannon, Region == region) ok = shapiro.test(sample$Shannon) valeur_shapiro = append(valeur_shapiro,ok$p.value) } data.frame(valeur_shapiro) # Pour nos deux régions (langue et mains) la distribution des valeurs de Shannon suit une distribution normale # Est ce que la variance est similaire ? # Nous pouvons facilement tester la vairance avec la fonction var.test var.test(Shannon ~ Region, shannon, alternative = "two.sided") t.test(Shannon ~ Region, data = shannon) Nos données respectent les postulats du test de t de Student, nous pouvons donc procéder à l’analyse. p = ggplot(shannon, aes(x = Region, y = Shannon, fill = Region)) + custom_theme() + geom_boxplot(alpha = 0.3, outlier.shape = NA) + geom_jitter(aes(color = Region), size = 3, ) + scale_fill_manual(values = c("cornflowerblue","palevioletred")) + scale_color_manual(values = c("cornflowerblue","palevioletred")) + labs(x = "Région", y = expression(paste("Indice de diversité de Shannon (", italic("H'"), ")"))) + scale_y_continuous(limits = c(0,4), expand = c(0,0)) + annotate(geom ="text", x = 1.5, y = 3.8, label = expression(paste(italic("p")," < .001")), size = 5) 5.3 Ordination (PCoA) # Transformer les données d'abondance relative avec l'indice d'Hellinger asv_hellinger = decostand(asv, method = "hellinger") # On combine les tableaux de donnnées dans un object phyloseq ps = phyloseq(otu_table(asv_hellinger, taxa_are_rows = TRUE), tax_table(as.matrix(taxa)), sample_data(meta)) # Calculer l'indice de distance de Bray-Curtis entre chaque échantillon dist = vegdist(t(asv_hellinger), method = "bray") # Générer l'ordination avec les fonction R de bases. PCOA = cmdscale(dist, eig=TRUE, add=TRUE) # On extrait les coordonnées des points position = PCOA$points # Changer le nom des colonnes colnames(position) = c("Axe.1", "Axe.2") # get percentage of variation explained by each axis percent_explained = 100*PCOA$eig/sum(PCOA$eig) # reduce number of digits (arrondir) reduced_percent = format(round(percent_explained[1:2], digits = 2), nsmall = 1, trim = TRUE) # Generate pretty labels for plot pretty_labs = c(glue("Axe 1 ({reduced_percent[1]}%)"), glue("Axe 2 ({reduced_percent[2]}%)")) # combine PCOA results with metadata df = merge(position, meta, by = 0) Générer le graphique plot = ggplot(df, aes(x=Axe.1, y=Axe.2, color=Region)) + custom_theme() + geom_point(size = 3) + labs(x = pretty_labs[1], y = pretty_labs[2], color = "Région") + scale_y_continuous(limits = c(-0.50,0.50), expand = c(0,0)) + scale_x_continuous(limits = c(-0.5,0.5), expand = c(0,0)) + geom_hline(yintercept = 0, linetype = "dashed", linewidth = 0.2) + geom_vline(xintercept = 0, linetype = "dashed", linewidth = 0.2) + scale_color_manual(values = c("cornflowerblue","palevioletred")) PERMANOVA et betadisper # Calculer la dispersion (distance au centroïde) dispr = vegan::betadisper(dist, group = meta$Region) boxplot(dispr, main = "", xlab = "") # ANOVA sur la dispersion (p = 0.001) permutest(dispr) # ANOVA sur la distance en fonction des région (p = 0.001) adns = adonis2(dist ~ meta$Region) 5.4 LEfSe # Définir le rang taxonomique d'intérêt taxa = "Genus" # On combine les tableaux de données dans un object phyloseq ps = phyloseq(otu_table(asv, taxa_are_rows=TRUE), tax_table(as.matrix(taxa)), sample_data(meta)) # Lancer l'analyse LEfSe out_lefse = run_lefse(ps, group = "Region",taxa_rank = taxa) # Extraire les résultats sous forme de tableau df_lefse = marker_table(out_lefse) %>% data.frame() # Keeping only the 9 taxa with the highest LDA score in each group number_of_taxa = 5 lda_out = list() i = 0 for (each_region in unique(df_lefse$enrich_group)){ i = i + 1 region = subset(df_lefse, enrich_group == each_region) top_taxa=head(region[order(region$ef_lda, decreasing= T),], n=number_of_taxa) lda_out[[i]] = top_taxa } df = do.call("rbind", lda_out) # Multiplier par -1 les score LDA pour la région sous la langue df$ef_lda = with(df, ifelse(enrich_group == "Sous la langue", -1 * df$ef_lda, 1*df$ef_lda)) Générer le graphique p = ggplot(df, aes(x = ef_lda, y = reorder(feature, ef_lda), fill = enrich_group)) + geom_bar(stat = "identity", width = 0.7, size = 0.5) + custom_theme() + facet_grid(rows = vars(enrich_group), scales = "free_y", space = "free_y") + theme(legend.position = "none", strip.text.y = element_text(angle =)) + scale_x_continuous(limits=c(-6,6), expand = c(0,0)) + labs(x = "Score LDA", y = "Genre") + scale_fill_manual(values=c("cornflowerblue","palevioletred")) + geom_vline(xintercept = 0, linewidth = 0.25, colour="black") 5.5 Abondances relatives # Defining the number of most abundant taxa to keep number_of_taxa = 3 ps_rel_abund=transform_sample_counts(ps, function(x) x/sum(x)) ps_glom = tax_glom(ps_rel_abund, taxrank = taxa_rank) melted_df = psmelt(ps_glom) # Create an empty list that we will populated with the unique taxa of each sample list_of_all_taxonomic_rank= list() i = 0 # Beginning of the for loop for (each_sample in unique(melted_df$Sample)){ i = i + 1 sample = subset(melted_df, Sample == each_sample) # Create a new dataframe from the iterator (sample). total_abundance = aggregate(sample$Abundance, by = list(taxa_rank = sample[[taxa_rank]]), FUN = sum) # Combine together the same taxa and sum the abundances top = head(total_abundance[order(total_abundance$x, decreasing= T),], n = number_of_taxa) # Sort by abundance and keep only the X number of taxa defined by variable number_of_taxa others_df = sample[!sample[[taxa_rank]] %in% top$taxa_rank,] # Extract in a new dataframe all taxa that are not present in the dataframe `top` others_list = pull(others_df, taxa_rank) # Create a list by pulling all the values from the column corresponding to the taxa_rank into a list sample[sample[[taxa_rank]]%in% others_list,][[taxa_rank]] = "Others" # In the dataframe `sample` rename all the taxa from the list `others_list` as `Others` list_of_all_taxonomic_rank[[i]] = sample #save this dataframe in our list } df = do.call("rbind",list_of_all_taxonomic_rank) # combine all the dataframe from the list into one dataframe unique_taxon = data.frame(unique(df[[taxa_rank]])) # create dataframe with the unique names of taxa name = colnames(unique_taxon) # extract the name of the column in order to rename the column with the following line names(unique_taxon)[names(unique_taxon) == name] = as.character(taxa_rank) # Rename the column to the taxa rank defined earlier # get the total number of unique most abundant taxa n = nrow(unique_taxon) # generate a set of X unique colors corresponding to the number of unique taxa palette = distinctColorPalette(n) unique_taxon[[taxa_rank]] = factor(unique_taxon[[taxa_rank]]) names(palette) = levels(unique_taxon[[taxa_rank]]) # assign gray to category "Others". The same nomenclature can be use to manually change certain colors. palette[["Others"]] = "#E1E1E1" # recreate palette with markdown to italicize name and remove the underscore after Unclassified all_names = data.frame(names(palette)) names_markdown = all_names %>% mutate(names.palette. = str_replace(names.palette., "(.*)","*\\\\1*"), # Adding asterisk at beginning and end of every taxa names.palette. = str_replace(names.palette., "\\\\*Unclassified_(.*)\\\\*","Unclassified *\\\\1*"), # Removing the asterisk for words that don't need to be italicize (Unclassified and Others) names.palette. = str_replace(names.palette., "\\\\*Others\\\\*", "Others")) list_names=as.vector(names_markdown$names.palette.) # Replace names of object names(palette) = c(list_names) # Making the same modification to the taxa name from the legend to the taxa names in the dataframe df[[taxa_rank]] = str_replace(df[[taxa_rank]], "(.*)","*\\\\1*") df[[taxa_rank]] = str_replace(df[[taxa_rank]], "\\\\*Unclassified_(.*)\\\\*","Unclassified *\\\\1*") df[[taxa_rank]] = str_replace(df[[taxa_rank]], "\\\\*Others\\\\*", "Others") # Ordering the legend in alphabetical order legend_raw = unique(df[[taxa_rank]]) #Extract legend as text ordered_legend = sort(legend_raw) # order alphabetically reordered_legend = fct_relevel(ordered_legend, "Others") # move "Others" to the beginning final_legend = levels(reordered_legend) # Extract the levels in a new object my_scale = scale_fill_manual(name = as.character(taxa_rank), breaks = paste(final_legend), values = palette, na.translate=FALSE, drop=TRUE, limits = force) p = ggplot(df, aes(x = Scientifque, weight = Abundance, fill = fct_reorder(.data[[taxa_rank]],Abundance,.desc=FALSE))) + # .data is very important to force the evaluation of the input variables (taxonomic_rank) geom_bar() + coord_flip() + labs(y ='Abondance relative', x="Scientifque") + theme_classic() + facet_wrap(~ Region, ncol=2, strip.position="top",scales='free_x') + theme(text = element_text(size = text_size), panel.spacing = unit(2, "lines"), plot.title = element_text(hjust =0.5, size=text_size), axis.title=element_text(size=text_size), axis.text.x = element_text(size=text_size, hjust=0.6), axis.text.y = element_text(size=text_size, vjust=0.5), legend.position = "bottom", legend.title = element_text(size = text_size), legend.text = element_markdown(size = 12), legend.key.size = unit(0.5, 'cm'), legend.margin = margin(), # pre-emptively set zero margins strip.background = element_blank(), strip.text = element_text(size = text_size)) + # remove facet_grid box background my_scale + # Load our color palette scale_y_continuous(limits = c(0, 1),breaks = seq(0, 1, by = 0.5), expand = c(0,0), labels = c("0", "0.5", "1")) + # Remove the white space # Adjusting the legend, notably the number of rows and position guides(fill = guide_legend(nrow = 8, title = "Genre", title.position = "top", title.hjust = 0.5, reverse=FALSE)) "],["404.html", "Page not found", " Page not found The page you requested cannot be found (perhaps it was moved or renamed). You may want to try searching to find the page's new location, or use the table of contents to find the page you are looking for. "]]