diff --git a/DE_varroa.Rmd b/DE_varroa.Rmd new file mode 100644 index 0000000..2c99ce1 --- /dev/null +++ b/DE_varroa.Rmd @@ -0,0 +1,419 @@ +--- +title: "Varroa microbiome" +author: "Sasha Mikheyev" +date: "1/23/2019" +output: html_document +--- + +```{r setup, include=FALSE} +knitr::opts_chunk$set(echo = TRUE, cache = T) +``` + +```{r library} +library(tidyverse) +library(sleuth) +``` + + +```{r sleuth} +bacteria <- read_tsv("mondet/gene_calls_nr.bacteria.tsv", col_names = c("contig", "protein_id", "taxname"), col_types = "-c---c-c") %>% filter(grepl("Bacteria", taxname)) %>% mutate(wolbachia = ifelse(grepl("Wolbachia", taxname), "yes", "no")) %>% separate(taxname, sep = ";\\s+", into = c("Bacteria",LETTERS[1:6]), extra = 'drop') %>% mutate(F = sub(";", "", F)) +samples <- tibble(id = dir(file.path("mondet/kallisto"))) %>% mutate(path = paste0("mondet/kallisto/", id)) %>% + left_join(read_tsv("mondet/treatments.tsv", col_types = cols())) %>% mutate(sample = id, condition = shortName) +``` + + +### Comparing males vs young + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("young", "males")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) + + +``` +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionyoung') +results_table = sleuth_results(so, 'conditionyoung') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +###Young vs Phoretic + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("young", "phoretic")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionyoung') +results_table = sleuth_results(so, 'conditionyoung') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Phoretics vs caged + +There were nio + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "phoreticCaged")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionphoreticCaged') +results_table = sleuth_results(so, 'conditionphoreticCaged') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` +### Phoretic vs arresting + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "arresting")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionphoretic') +results_table = sleuth_results(so, 'conditionphoretic') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` + +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` +### Аrresting vs. pre-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "prelaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionprelaying') +results_table = sleuth_results(so, 'conditionprelaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` +### Pre-laying vs. Laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "laying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionphoretic') +results_table = sleuth_results(so, 'conditionphoretic') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` +### Laying vs. post-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("laying", "postlaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionpostlaying') +results_table = sleuth_results(so, 'conditionpostlaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Post-laying vs. Emerging + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("postlaying", "emerging")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) + +sleuth_significant %>% left_join(bacteria, by = c("target_id" = "contig")) %>% dplyr::filter(!is.na(A)) %>% ggplot(aes(A))+geom_histogram(stat = "count") + scale_y_log10() + coord_flip() +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionpostlaying') +results_table = sleuth_results(so, 'conditionpostlaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Post-laying vs. non-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("postlaying", "nonlaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionpostlaying') +results_table = sleuth_results(so, 'conditionpostlaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Laying vs. non-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("laying", "nonlaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionnonlaying') +results_table = sleuth_results(so, 'conditionnonlaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Emerging vs. Phoretic + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("emerging", "phoretic")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionphoretic') +results_table = sleuth_results(so, 'conditionphoretic') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Emerging vs. Young + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("emerging", "young")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionyoung') +results_table = sleuth_results(so, 'conditionyoung') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` \ No newline at end of file diff --git a/DE_varroa_microbiome.csv b/DE_varroa_microbiome.csv new file mode 100644 index 0000000..83e31e2 --- /dev/null +++ b/DE_varroa_microbiome.csv @@ -0,0 +1,13 @@ +Stages,M_DE,M_up-regulated,M_down-regulated,V_DE,V_up-regulated,V_down-regulated,V2_DE,V2_up-regulated,V2_down-regulated +males vs young,5397,1596,3801,19220,8328,10892,8069,2641,5428 +Young vs Phoretic,2852,1570,1282,12769,6846,5923,3374,1763,1611 +Phoretics vs caged,1,0,1,0,0,0,23,13,10 +Phoretic vs arresting,68,19,49,3514,1288,2226,100,36,64 +Arrest vs pre-laying,2866,1527,1339,9458,4769,4689,3715,1890,1825 +Pre-laying vs. Laying,4247,1951,2296,2683,1795,888,5666,2654,3012 +Laying vs. post-laying,2803,978,1105,9456,4077,5379,2851,1264,1587 +Post-laying vs. Emerging,619,322,297,3536,2135,1401,619,322,297 +Post-laying vs. Non-laying,5,3,2,2,2,0,4,1,3 +Laying vs Non-laying,2203,1903,1110,13065,5738,7327,2203,1093,1110 +Emerging vs Phoretic,571,219,352,1492,889,603,571,219,352 +Emerging vs Young,1870,1011,859,8947,4465,4482,1870,1011,859 \ No newline at end of file diff --git a/DE_varroa_mite.Rmd b/DE_varroa_mite.Rmd new file mode 100644 index 0000000..513602c --- /dev/null +++ b/DE_varroa_mite.Rmd @@ -0,0 +1,417 @@ +--- +title: "Varroa" +author: "Lazzat Aibekova" +date: "6/04/2019" +output: html_document +--- + +```{r setup, include=FALSE} +knitr::opts_chunk$set(echo = TRUE, cache = T) +``` + +```{r library} +library(tidyverse) +library(sleuth) +``` + + +```{r sleuth} +samples <- tibble(id = dir(file.path("mondet/kallisto_mite"))) %>% mutate(path = paste0("mondet/kallisto_mite/", id)) %>% + left_join(read_tsv("mondet/treatments.tsv", col_types = cols())) %>% mutate(sample = id, condition = shortName) +``` + + +### Comparing males vs young + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("young", "males")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) + + +``` +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionyoung') +results_table = sleuth_results(so, 'conditionyoung') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +###Young vs Phoretic + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("young", "phoretic")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionyoung') +results_table = sleuth_results(so, 'conditionyoung') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Phoretics vs caged + +There were nio + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "phoreticCaged")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionphoreticCaged') +results_table = sleuth_results(so, 'conditionphoreticCaged') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` +### Phoretic vs arresting + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "arresting")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionphoretic') +results_table = sleuth_results(so, 'conditionphoretic') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` + +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` +### Аrresting vs. pre-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "prelaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionprelaying') +results_table = sleuth_results(so, 'conditionprelaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` +### Pre-laying vs. Laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "laying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionphoretic') +results_table = sleuth_results(so, 'conditionphoretic') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` +### Laying vs. post-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("laying", "postlaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionpostlaying') +results_table = sleuth_results(so, 'conditionpostlaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Post-laying vs. Emerging + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("postlaying", "emerging")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) + +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionpostlaying') +results_table = sleuth_results(so, 'conditionpostlaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Post-laying vs. non-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("postlaying", "nonlaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionpostlaying') +results_table = sleuth_results(so, 'conditionpostlaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Laying vs. non-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("laying", "nonlaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionnonlaying') +results_table = sleuth_results(so, 'conditionnonlaying') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Emerging vs. Phoretic + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("emerging", "phoretic")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionphoretic') +results_table = sleuth_results(so, 'conditionphoretic') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` + +### Emerging vs. Young + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("emerging", "young")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +```{r} +models(so) +``` +Number of Differentially expressed genes (TRUE): +```{r} +so = sleuth_wt(so, 'conditionyoung') +results_table = sleuth_results(so, 'conditionyoung') +Number_DE_genes = table(results_table$qval <= 0.05) +Number_DE_genes +``` +Number of up-regulated (TRUE) and down-regulated (FALSE) genes: +```{r} +results_table +DE_genes = subset(results_table, qval <= 0.05) +DE_genes +up_down_regulated = table(DE_genes$b > 0) +up_down_regulated + +``` \ No newline at end of file diff --git a/Plots.Rmd b/Plots.Rmd new file mode 100644 index 0000000..8f292bc --- /dev/null +++ b/Plots.Rmd @@ -0,0 +1,39 @@ +--- +title: "R Notebook" +output: html_notebook +--- +```{r} +data = read.csv("/Users/lazzataibekova/Dropbox (OIST)/varroa microbiome/DE_varroa_microbiome.csv") +data + +``` + +```{r} +plot(data$V_DE, data$M_DE, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Total number of DEG") +cor(data$V_DE, data$M_DE) +``` + +```{r} +plot(data$V_up.regulated, data$M_up.regulated, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Number of up-regulated genes") +cor(data$V_up.regulated, data$M_up.regulated) +``` + +```{r} +plot(data$V_down.regulated, data$M_down.regulated, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Number of down-regulated genes") +cor(data$V_down.regulated, data$M_down.regulated) +``` + +```{r} +plot(data$V2_DE, data$M_DE, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Total number of DEG (re-mapped)") +cor(data$V2_DE, data$M_DE) +``` +```{r} +plot(data$V2_up.regulated, data$M_up.regulated, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Number of up-regulated genes (re-mapped)") +cor(data$V2_up.regulated, data$M_up.regulated) +``` + +```{r} +plot(data$V2_down.regulated, data$M_down.regulated, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Number of down-regulated genes (re-mapped)") +cor(data$V2_down.regulated, data$M_down.regulated) +``` + diff --git a/Plots.nb.html b/Plots.nb.html new file mode 100644 index 0000000..6fbf60a --- /dev/null +++ b/Plots.nb.html @@ -0,0 +1,399 @@ + + + + + + + + + + + + + +R Notebook + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +
+ + + + + + + + + + + +
data = read.csv("/Users/lazzataibekova/Dropbox (OIST)/varroa microbiome/DE_varroa_microbiome.csv")
+data
+ + +
+ +
+ + + + + + +
plot(data$V_DE, data$M_DE, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Total number of DEG")
+ + +

+ + +
cor(data$V_DE, data$M_DE)
+ + +
[1] 0.7540207
+ + + + + + +
plot(data$V_up.regulated, data$M_up.regulated, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Number of up-regulated genes")
+ + +

+ + +
cor(data$V_up.regulated, data$M_up.regulated)
+ + +
[1] 0.7598868
+ + + + + + +
plot(data$V_down.regulated, data$M_down.regulated, xlab = "Varroa genes", ylab = "Microbiome genes", main = "Number of down-regulated genes")
+ + +

+ + +
cor(data$V_down.regulated, data$M_down.regulated)
+ + +
[1] 0.7331989
+ + + + + + +
plot(data$V2_DE, data$M_DE,  xlab = "Varroa genes", ylab = "Microbiome genes", main = "Total number of DEG (re-mapped)")
+ + +

+ + +
cor(data$V2_DE, data$M_DE)
+ + +
[1] 0.9840727
+ + + + +
plot(data$V2_up.regulated, data$M_up.regulated,  xlab = "Varroa genes", ylab = "Microbiome genes", main = "Number of up-regulated genes (re-mapped)")
+ + +

+ + +
cor(data$V2_up.regulated, data$M_up.regulated)
+ + +
[1] 0.898323
+ + + + + + +
plot(data$V2_down.regulated, data$M_down.regulated,  xlab = "Varroa genes", ylab = "Microbiome genes", main = "Number of down-regulated genes (re-mapped)")
+ + +

+ + +
cor(data$V2_down.regulated, data$M_down.regulated)
+ + +
[1] 0.9947204
+ + + + + +
LS0tCnRpdGxlOiAiUiBOb3RlYm9vayIKb3V0cHV0OiBodG1sX25vdGVib29rCi0tLQpgYGB7cn0KZGF0YSA9IHJlYWQuY3N2KCIvVXNlcnMvbGF6emF0YWliZWtvdmEvRHJvcGJveCAoT0lTVCkvdmFycm9hIG1pY3JvYmlvbWUvREVfdmFycm9hX21pY3JvYmlvbWUuY3N2IikKZGF0YQoKYGBgCgpgYGB7cn0KcGxvdChkYXRhJFZfREUsIGRhdGEkTV9ERSwgeGxhYiA9ICJWYXJyb2EgZ2VuZXMiLCB5bGFiID0gIk1pY3JvYmlvbWUgZ2VuZXMiLCBtYWluID0gIlRvdGFsIG51bWJlciBvZiBERUciKQpjb3IoZGF0YSRWX0RFLCBkYXRhJE1fREUpCmBgYAoKYGBge3J9CnBsb3QoZGF0YSRWX3VwLnJlZ3VsYXRlZCwgZGF0YSRNX3VwLnJlZ3VsYXRlZCwgeGxhYiA9ICJWYXJyb2EgZ2VuZXMiLCB5bGFiID0gIk1pY3JvYmlvbWUgZ2VuZXMiLCBtYWluID0gIk51bWJlciBvZiB1cC1yZWd1bGF0ZWQgZ2VuZXMiKQpjb3IoZGF0YSRWX3VwLnJlZ3VsYXRlZCwgZGF0YSRNX3VwLnJlZ3VsYXRlZCkKYGBgCgpgYGB7cn0KcGxvdChkYXRhJFZfZG93bi5yZWd1bGF0ZWQsIGRhdGEkTV9kb3duLnJlZ3VsYXRlZCwgeGxhYiA9ICJWYXJyb2EgZ2VuZXMiLCB5bGFiID0gIk1pY3JvYmlvbWUgZ2VuZXMiLCBtYWluID0gIk51bWJlciBvZiBkb3duLXJlZ3VsYXRlZCBnZW5lcyIpCmNvcihkYXRhJFZfZG93bi5yZWd1bGF0ZWQsIGRhdGEkTV9kb3duLnJlZ3VsYXRlZCkKYGBgCgpgYGB7cn0KcGxvdChkYXRhJFYyX0RFLCBkYXRhJE1fREUsICB4bGFiID0gIlZhcnJvYSBnZW5lcyIsIHlsYWIgPSAiTWljcm9iaW9tZSBnZW5lcyIsIG1haW4gPSAiVG90YWwgbnVtYmVyIG9mIERFRyAocmUtbWFwcGVkKSIpCmNvcihkYXRhJFYyX0RFLCBkYXRhJE1fREUpCmBgYApgYGB7cn0KcGxvdChkYXRhJFYyX3VwLnJlZ3VsYXRlZCwgZGF0YSRNX3VwLnJlZ3VsYXRlZCwgIHhsYWIgPSAiVmFycm9hIGdlbmVzIiwgeWxhYiA9ICJNaWNyb2Jpb21lIGdlbmVzIiwgbWFpbiA9ICJOdW1iZXIgb2YgdXAtcmVndWxhdGVkIGdlbmVzIChyZS1tYXBwZWQpIikKY29yKGRhdGEkVjJfdXAucmVndWxhdGVkLCBkYXRhJE1fdXAucmVndWxhdGVkKQpgYGAKCmBgYHtyfQpwbG90KGRhdGEkVjJfZG93bi5yZWd1bGF0ZWQsIGRhdGEkTV9kb3duLnJlZ3VsYXRlZCwgIHhsYWIgPSAiVmFycm9hIGdlbmVzIiwgeWxhYiA9ICJNaWNyb2Jpb21lIGdlbmVzIiwgbWFpbiA9ICJOdW1iZXIgb2YgZG93bi1yZWd1bGF0ZWQgZ2VuZXMgKHJlLW1hcHBlZCkiKQpjb3IoZGF0YSRWMl9kb3duLnJlZ3VsYXRlZCwgZGF0YSRNX2Rvd24ucmVndWxhdGVkKQpgYGAKCg==
+ + + +
+ + + + + + + + diff --git a/mondet.Rmd b/mondet.Rmd new file mode 100644 index 0000000..0105273 --- /dev/null +++ b/mondet.Rmd @@ -0,0 +1,136 @@ +--- +title: "Varroa microbiome" +author: "Sasha Mikheyev" +date: "1/23/2019" +output: html_document +--- + +```{r setup, include=FALSE} +knitr::opts_chunk$set(echo = TRUE, cache = T) +``` + +```{r library} +library(tidyverse) +library(sleuth) +``` + + +```{r sleuth} +bacteria <- read_tsv("mondet/gene_calls_nr.bacteria.tsv", col_names = c("contig", "protein_id", "taxname"), col_types = "-c---c-c") %>% filter(grepl("Bacteria", taxname)) %>% mutate(wolbachia = ifelse(grepl("Wolbachia", taxname), "yes", "no")) %>% separate(taxname, sep = ";\\s+", into = c("Bacteria",LETTERS[1:6]), extra = 'drop') %>% mutate(F = sub(";", "", F)) +samples <- tibble(id = dir(file.path("mondet/kallisto"))) %>% mutate(path = paste0("mondet/kallisto/", id)) %>% + left_join(read_tsv("mondet/treatments.tsv", col_types = cols())) %>% mutate(sample = id, condition = shortName) +``` + +### Comparing males vs young + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("young", "males")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) + + +``` + + +### Phoretics vs caged + +There were nio + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "phoreticCaged")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + + +### Phoretic vs arresting + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "arresting")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +### Аrresting vs. pre-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "prelaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + + +### Pre-laying vs. Laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("phoretic", "laying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +### Laying vs. post-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("laying", "postlaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` + +### Post-laying vs. Emerging + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("postlaying", "emerging")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) + +sleuth_significant %>% left_join(bacteria, by = c("target_id" = "contig")) %>% dplyr::filter(!is.na(A)) %>% ggplot(aes(A))+geom_histogram(stat = "count") + scale_y_log10() + coord_flip() +``` + + +### Post-laying vs. non-laying + +```{r} +so <- sleuth_prep(dplyr::filter(samples, shortName %in% c("postlaying", "nonlaying")), extra_bootstrap_summary = TRUE) +so <- sleuth_fit(so, ~condition, 'full') +so <- sleuth_fit(so, ~1, 'reduced') +so <- sleuth_lrt(so, 'reduced', 'full') +sleuth_table <- sleuth_results(so, 'reduced:full', 'lrt', show_all = FALSE) +sleuth_significant <- dplyr::filter(sleuth_table, qval <= 0.05) +dim(sleuth_significant) +head(sleuth_significant) +``` \ No newline at end of file