Load Packages

library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.2.1     ✔ readr     2.2.0
## ✔ forcats   1.0.1     ✔ stringr   1.6.0
## ✔ ggplot2   4.0.3     ✔ tibble    3.3.1
## ✔ lubridate 1.9.5     ✔ tidyr     1.3.2
## ✔ purrr     1.2.2     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

Import and Investigate Files

conditionData <- read.csv(
  file = "conditions_annotation.csv",
  header = TRUE,
  stringsAsFactors = FALSE,
  check.names = FALSE
)

expressionData <- read.csv(
  file = "SC_expression.csv",
  header = TRUE,
  stringsAsFactors = FALSE,
  check.names = FALSE
)

labelData <- read.csv(
  file = "mergedLabels.csv",
  header = TRUE,
  stringsAsFactors = FALSE
)

names(conditionData)
names(expressionData)
names(labelData)

unique(conditionData$primary)

Rename expressionData column 1

names(expressionData)[1] <- "gene"

Choose condition and filter conditionData

unique(conditionData$primary)

swr1 <- conditionData[
  grepl("swr1", conditionData$primary),
]
names(swr1)

Extract treatments and filter

treatments <- swr1$ID
treatments

filter <- expressionData %>%
  select(gene, all_of(treatments))

Merge filtered expressionData with labelData

merged <- inner_join(labelData,
                     filter,
                     by = "gene")

Pivot Longer

mergedLong <- merged %>%
  pivot_longer(
    cols = all_of(treatments),
    names_to = "treatment",
    values_to = "count"
  )

Tibble of mean and median

treatmentSum <- mergedLong %>%
  group_by(treatment) %>%
  summarise(
    mean = mean(count, na.rm = TRUE),
    median = median(count, na.rm = TRUE),
    n = n()
  )

treatmentSum
## # A tibble: 2 × 4
##   treatment  mean median     n
##   <chr>     <dbl>  <dbl> <int>
## 1 AFNCCA     37.8  10.1   6011
## 2 AFNCCR     17.8   4.77  6011

Filter Outliers

avgmeans <- mean(mergedLong$count, na.rm = TRUE)
avgmeans
## [1] 27.76336
cutoff <- 2 * avgmeans
cutoff
## [1] 55.52671
plotData <- mergedLong %>%
  filter(count <= cutoff)

Plot

ggplot(plotData,
       aes(x = treatment,
           y = count)) +
  geom_violin() +
  labs(
    title = "Distribution: Gene Expression by Treatment",
    x = "Treatment",
    y = "Expression Count",
  )