---
title: "simulate_reads_summary"
author: "bhaas"
date: '2024-02-01'
output: github_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
library(tidyverse)
```


Here we build a summary figure for the pacbio HiFi and ONT R10.4.1 chemistry simulated reads (via pbsim3).

See the subdirectories for the individual analyses:
    1b.1.PacBio_simulations 
       and
    1b.2.ONT_simulations


```{r}

prog_ordering = c('ctat-LR-fusion', 'JAFFAL', 'LongGF', 'fusionseeker', 'pbfusion')


```



# Examine peak F1 score

```{r}

pacbio_F1_summary_data = read.table("1b.1.PacBio_simulations/PacBio_sim.mean_F1.combined_results.tsv", header=T, sep="\t")

pacbio_F1_summary_data = pacbio_F1_summary_data %>% mutate(seqTech = "PacBio")

pacbio_F1_summary_data$prog = factor(pacbio_F1_summary_data$prog, levels=prog_ordering)

pacbio_F1_summary_data %>% head()
```



```{r}

ONT_F1_summary_data = read.table("1b.2.ONT_simulations/ONT_sim.mean_F1.combined_results.tsv", header=T, sep="\t")

ONT_F1_summary_data = ONT_F1_summary_data %>% mutate(seqTech = "ONT")

ONT_F1_summary_data %>% head()

```


```{r}

combined_F1_data = bind_rows(pacbio_F1_summary_data, ONT_F1_summary_data)


combined_F1_data %>% head()
```


```{r}

combined_F1_data$prog = factor(combined_F1_data$prog, levels = c('ctat-LR-fusion',
                                                           'JAFFAL',
                                                           'LongGF',
                                                           'fusionseeker',
                                                           'pbfusion') )


```




```{r}


combined_F1_data %>%
    mutate(analysisType = factor(analysisType, 
                            levels=c('strict', 'allow_reverse', 'Exact Brkpt', 'Fuzzy Brkpt') )) %>%
    mutate(seqTech = factor(seqTech, levels=c("PacBio", "ONT"))) %>%
     ggplot() + theme_bw() +
     geom_jitter(aes(x=analysisType, y=mean_F1, color=prog, shape=prog), width=0.2, height=0, size=rel(3)) +
     facet_grid(vars(seqTech), vars(coverage_level)) +
     theme(axis.text.x = element_text(angle = 90, hjust = 1)) 


```





```{r}

combined_F1_paperfig = combined_F1_data %>%
    filter(coverage_level=="cov5") %>%
    mutate(analysisType = factor(analysisType, 
                            levels=c('strict', 'allow_reverse', 'Exact Brkpt', 'Fuzzy Brkpt') )) %>%
    mutate(seqTech = factor(seqTech, levels=c("PacBio", "ONT"))) %>%
     ggplot() + theme_bw() +
     geom_jitter(aes(x=analysisType, y=mean_F1, color=prog, shape=prog), width=0.2, height=0, size=rel(3)) +
     facet_wrap(~seqTech) +
     theme(axis.text.x = element_text(angle = 90, hjust = 1)) 


combined_F1_paperfig




```

```{r}

ggsave(combined_F1_paperfig, filename="pbsim3_bmark.F1.paperfig.svg", width=8, height=5)


```


# Examine PR-AUC


```{r}

pacbio_PR_AUC_summary_data = read.table("1b.1.PacBio_simulations/PacBio_sim.mean_PR_AUC.combined_results.tsv", header=T, sep="\t")

pacbio_PR_AUC_summary_data = pacbio_PR_AUC_summary_data %>% mutate(seqTech = "PacBio")


pacbio_PR_AUC_summary_data %>% head()
```



```{r}

ONT_PR_AUC_summary_data = read.table("1b.2.ONT_simulations/ONT_sim.mean_PR_AUC.combined_results.tsv", header=T, sep="\t")

ONT_PR_AUC_summary_data = ONT_PR_AUC_summary_data %>% mutate(seqTech = "ONT")

ONT_PR_AUC_summary_data %>% head()

```


```{r}

combined_PR_AUC_data = bind_rows(pacbio_PR_AUC_summary_data, ONT_PR_AUC_summary_data)


combined_PR_AUC_data %>% head()
```


```{r}

combined_PR_AUC_data$prog = factor(combined_PR_AUC_data$prog, levels = prog_ordering )

```




```{r}


combined_PR_AUC_data %>%
    mutate(analysisType = factor(analysisType, 
                            levels=c('strict', 'allow_reverse', 'Exact Brkpt', 'Fuzzy Brkpt') )) %>%
    mutate(seqTech = factor(seqTech, levels=c("PacBio", "ONT"))) %>%
     ggplot() + theme_bw() +
     geom_jitter(aes(x=analysisType, y=mean_AUC, color=prog, shape=prog), width=0.2, height=0, size=rel(3)) +
     facet_grid(vars(seqTech), vars(coverage_level)) +
     theme(axis.text.x = element_text(angle = 90, hjust = 1)) 


```





```{r}

combined_PR_AUC_paperfig = combined_PR_AUC_data %>%
    filter(coverage_level=="cov5") %>%
    mutate(analysisType = factor(analysisType, 
                            levels=c('strict', 'allow_reverse', 'Exact Brkpt', 'Fuzzy Brkpt') )) %>%
    mutate(seqTech = factor(seqTech, levels=c("PacBio", "ONT"))) %>%
     ggplot() + theme_bw() +
     geom_jitter(aes(x=analysisType, y=mean_AUC, color=prog, shape=prog), width=0.2, height=0, size=rel(3)) +
     facet_wrap(~seqTech) +
     theme(axis.text.x = element_text(angle = 90, hjust = 1)) 


combined_PR_AUC_paperfig




```

```{r}

ggsave(combined_PR_AUC_paperfig, filename="pbsim3_bmark.AUC.paperfig.svg", width=8, height=5)


```

