# Differential expression and abundance between conditions

```{r multiSplComp_setup, include=FALSE, echo=FALSE}
# First, set some variables:
require(knitr)

opts_chunk$set(error=FALSE, message=FALSE, warning=FALSE)
opts_chunk$set(echo = TRUE)
opts_chunk$set(eval = TRUE) 
options(stringsAsFactors = FALSE)
opts_chunk$set(fig.width=7, fig.height=7)
set.seed(123) # for reproducibility
```

```{r}
splSetToGet <- "PBMMC,ETV6-RUNX1"
splSetVec <- unlist(strsplit(splSetToGet, ","))
splSetToGet2 <- gsub(",", "_", splSetToGet)
nbPcToComp <- 50
figSize <- 7
```

```{r, message=FALSE, warning=FALSE}
library(scater)
library(scran)
library(batchelor)
library(edgeR)
library(tidyverse)
library(patchwork)
library(DT)
fontsize <- theme(axis.text=element_text(size=12), axis.title=element_text(size=16))
```

Source: [Multi-sample comparisons](https://osca.bioconductor.org/multi-sample-comparisons.html) of the OSCA book.

## Exercise 2 - differential abundance

Imagine ETV6-RUNX1_4 failed, leaving you with three ETV6-RUNX1 replicates ... but all else remains as above, including the clusters identified.

Identify clusters whose abundance differ between conditions.

## Setting up the data

Load the SCE object (with 1200 cells per sample):

```{r}
# Read object in:
merged <- readRDS("../Robjects/caron_sce_nz_postDeconv_1p2kcps_dsi_PBMMC_ETV6-RUNX1_merged.Rds")
# also get raw counts that were written to a separate file
# (to help file sharing)
merged_counts <- readRDS("../Robjects/caron_sce_nz_postDeconv_1p2kcps_dsi_PBMMC_ETV6-RUNX1_merged_counts.Rds")
# put raw counts back:
counts(merged) <- merged_counts
# tidy:
rm(merged_counts)
```

```{r}
bcToKeep <- colData(merged) %>%
  data.frame() %>%
  rownames_to_column("Barcode") %>%
  filter(!SampleName == "ETV6-RUNX1_4") %>%
  pull(Barcode)
indToKeep <- which(colnames(merged) %in% bcToKeep)
merged <- merged[,indToKeep]

merged$SampleName <- factor(merged$SampleName)
```

A brief inspection of the results shows clusters contain varying contributions from samples:

```{r}
colLabels(merged) <- merged$clusters.mnn
tab <- table(colLabels(merged), merged$SampleName)
tab
```

Count cells assigned to each cluster:

```{r}


```


Prepare edgeR object

```{r}
# think 'DGEList'





```


Filter out low-abundance labels

```{r}
# think 'filterByExpr'


```

Make design matrix

```{r}
design <- model.matrix(~factor(SampleGroup), y.ab$samples)
```

Estimate the NB dispersion

```{r}
# think 'estimateDisp'

```

Estimate the QL dispersion

```{r}
# think 'glmQLFit'



```

Test for differences in abundance between sample groups using `glmQLFTest`

```{r}


```

Get summary

```{r}
# think 'decideTests'

```

Display whole table:

```{r}
# think 'topTags'

```


