---
title: "02-kallisto"
output: html_document
---

```{bash}
/home/shared/kallisto/kallisto \
index -i /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/genome/Panopea-generosa-v1.0.index \
/home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/genome/Panopea-generosa-v1.0.fa
```
```{bash}
/home/shared/kallisto/kallisto \
quant
```




```{bash}
find /home/shared/8TB_HDD_01/sr320/ncbi/*_1.fastq \
| xargs basename -s _1.fastq | xargs -I{} /home/shared/kallisto/kallisto \
quant -i /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/data/transcriptome_v5.idx \
-o /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/{} \
-t 40 \
/home/shared/8TB_HDD_01/sr320/ncbi/{}_1.fastq /home/shared/8TB_HDD_01/sr320/ncbi/{}_2.fastq \
2>&1 | tee /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/{}.out
```

```{bash}
find /home/shared/8TB_HDD_01/sr320/ncbi/*_1.fastq \
| xargs basename -s _1.fastq | xargs -I{} /home/shared/kallisto/kallisto \
quant -i /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/data/transcriptome_v5.idx \
-o /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/b_{} \
-t 40 \
-bias \
/home/shared/8TB_HDD_01/sr320/ncbi/{}_1.fastq \
/home/shared/8TB_HDD_01/sr320/ncbi/{}_2.fastq
```

```{bash}
perl /home/shared/trinityrnaseq-v2.12.0/util/abundance_estimates_to_matrix.pl \
--est_method kallisto \
    --gene_trans_map none \
    --out_prefix /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/kallisto-0206 \
    --name_sample_by_basedir \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12207404/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12207405/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12207406/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12207407/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12212519/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12218868/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12218869/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12218870/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12226692/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12226693/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12227929/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12227930/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR12227931/abundance.tsv \
     /home/shared/8TB_HDD_01/sr320/github/nb-2022/P_generosa/analyses/SRR8788211/abundance.tsv 
```     

```{r}
library(tidyverse)
```


```{r}
du <- read.delim("../analyses/kallisto-0206.isoform.counts.matrix", header = TRUE, sep = '\t')
```

```{r}
ggplot(du, aes(x = SRR12207404)) +
  geom_histogram() +
  scale_y_log10() +
  xlim(0, 100000) 
  #ylim(0, 100)
```

```{r}
ggplot(du, aes(x = SRR12207405)) +
  geom_histogram() +
  scale_y_log10() +
  xlim(0, 100000) 
  #ylim(0, 100)
```

```{r}
ggplot(du, aes(x = SRR12207406)) +
  geom_histogram(bins = 100) +
  scale_y_log10() +
  xlim(0, 100000) 
  #ylim(0, 100)
```

```{r}
ggplot(du, aes(x = SRR12207407)) +
  geom_histogram(bins = 100) +
  scale_y_log10() +
  xlim(0, 100000) 
  #ylim(0, 100)
```



# plot them all

```{r}
dul <- pivot_longer(du, cols = starts_with("SRR"),
                       names_to = "library") 

```



```{r}
ggplot(data = dul, mapping = aes(x = value)) +
  geom_boxplot() +
  facet_wrap(~library)

```

```{r}
ggplot(dul, aes(library, value)) +
  geom_boxplot(color = "blue") +
  geom_jitter(color="black", size=0.4, alpha=0.9) 
  #theme(axis.text.x = element_text(angle = 90, vjust = 0.5))
```

```{r}
ggplot(dul, aes(library, value)) +
  geom_boxplot(color = "blue") +
  scale_y_log10() 
  #theme(axis.text.x = element_text(angle = 90, vjust = 0.5))
```




```{r}
ggplot(dul2, aes(x = value)) +
  geom_histogram(bins = 1000) +
  scale_y_log10() +
  facet_wrap(~library) +
  xlim(0, 10000) 
  #ylim(0, 100)
```
```{r}
ggplot(dul, aes(x=value, y=X)) + 
  geom_boxplot(aes(col=library)) + 
  ggtitle("Length by treatment (and preservation method)") + 
  theme_minimal() +
  geom_jitter(width=0.2, aes(col=library)) 
  #facet_wrap(~ethanol)
```

```{r}
library(hrbrthemes)
library(viridis)
```




```{r}
dul %>%
  ggplot( aes(x=library, y=value, fill=library)) +
    geom_boxplot() +
   # scale_fill_viridis(discrete = TRUE, alpha=0.6) +
    geom_jitter(color="black", size=0.4, alpha=0.9) +
   # theme_ipsum() +
    theme(
      legend.position="none",
      plot.title = element_text(size=11)
    ) +
    ggtitle("A boxplot with jitter") +
    xlab("")
```

```{r}
dul %>%
  ggplot( aes(x=library, y=value, fill=library)) +
    geom_violin() +
   # scale_fill_viridis(discrete = TRUE, alpha=0.6, option="A") +
    #theme_ipsum() +
    theme(
      legend.position="none",
      plot.title = element_text(size=11)
    ) +
    ggtitle("Violin chart") +
    ylim(0, 100) +
     scale_y_log10() +
    xlab("")

```
```{r}
library(tidyverse)

# Diamonds dataset is provided by R natively
#head(diamonds)
 
# basic example
dul %>%
  ggplot(aes(x = value, y = library)) +
  geom_density_ridges() +
  theme_ridges() +
  theme(legend.position = "none")
```

y=Count, fill=factor(Count))


```{r}
ggplot(dul, aes(y=value)) + 
  geom_bar(aes(fill=factor(value))) +
  facet_wrap(~library)
```

dat$col <- cut(dat$y,
               breaks = c(-Inf, 3000, Inf),
               labels = c("<=3000", ">3000"))

ggplot(dat, aes(x = x, y = y, color = col)) +
  geom_point()

```{r}
library(tidyverse)
```


```{r}
dul$col <- cut(dul$value, 
                breaks = c(-Inf, 0, 5, 10, 1000, Inf),
                labels = c("zero", "little", "lm", "1k", "high"))
```

```{r}
dut <- mutate(du, totexp = rowMeans(select(du, starts_with("SRR")), na.rm = TRUE)) %>%
  filter(totexp > 5)
```

```{r}
dul2 <- pivot_longer(dut, cols = starts_with("SRR"),
                       names_to = "library") 

```

mutate(df, IVMean = rowMeans(select(df, starts_with("IV")), na.rm = TRUE))


```{r}
dul2$col <- cut(dul2$totexp, 
                breaks = c(0, 5, 10, 1000, Inf),
                labels = c("little", "lm", "1k", "high"))
```

```{r}
ggplot(dul2, aes(factor(value), fill = factor(col))) + geom_bar()

```


```{r}
ggplot(dul, aes(factor(library), fill = factor(col))) + geom_bar()

```

```{r}
ggplot(dul2, aes(value, fill = factor(col))) + geom_bar() +
  facet_wrap(~library)

```

ggplot(acacia, aes(x = CIRC, y = HEIGHT)) +
  geom_point()




