---
title: "Population Analysis by Lineage - Non-Bottleneck Lineages related to Fig 1"
author: "Danielle Ireland"
date: "2025-06-24"
output: html_document
---

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


## Load libraries

```{r cars}
library(tidyverse)
library(dplyr)
library(readxl)
library(ggplot2)

set.seed(228)
```


## Read in non-bottleneck division data 

```{r }
div_NB_original<-as.data.frame(read_excel("Dataset_S1.xlsx", sheet="Dataset S1A"))
#select only the columns needed for further analysis
div_NB<-div_NB_original[c(10,11,12,14)]
#make a new column listing if this is a death (0=no, this is division). to be used to count the total number of worms in a generation
div_NB$Death=0
```

## Read in non-bottleneck DEATH data 
```{r}
dth_NB_original<-as.data.frame(read_excel("Dataset_S1.xlsx", sheet="Dataset S1B"))
#select only the columns needed for further analysis
dth_NB<-dth_NB_original[c(10:13)]
#make a new column listing if this is a death (1=YES). to be used to count the total number of worms in a generation
dth_NB$Death=1
```


#combine division and death data into one dataframe
```{r}
divdth <-rbind(div_NB, dth_NB)
divdth$FamilyName<-as.factor(divdth$FamilyName)
divdth$Type<-as.factor(ifelse(divdth$H.Gen!=0, 'H', 'T'))
#group by head and tail offspring
divdth_grouped<-split(divdth,list(divdth$FamilyName,divdth$Type))
```

#Determine death rate per lineage for heads and tail offspring separately (Fig 1E)
```{r}
#functions to count the number of alive and dead
dthNum<-function(x){
  colSums(x[5])}
AliveNum<-function(x){
  colSums(x[5]==0)}

NumDead<-sapply(divdth_grouped,dthNum)
NumAlive<-sapply(divdth_grouped,AliveNum)
DeathCount=as.data.frame(t(rbind(NumDead,NumAlive)))
DeathCount$FamilyName <- as.factor(c('F', 'H', 'P', 'S','F', 'H', 'P', 'S'))
DeathCount$Type <- as.factor(c('Head', 'Head', 'Head', 'Head', 'Tail', 'Tail', 'Tail', 'Tail'))
DeathCount$Percent<-100*DeathCount$NumDead/(DeathCount$NumDead+DeathCount$NumAlive)
DeathCount$Total<-DeathCount$NumDead+DeathCount$NumAlive

Fig1E<-ggplot(DeathCount, aes(x=(reorder(FamilyName, Percent)), y=Percent, fill=Type))+
  geom_bar(position='dodge', stat="identity", color='black')+
   scale_fill_grey()+
  xlab("Lineage")+
  ylab('% Death')+
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
        panel.background = element_blank(), axis.line = element_line(colour = "black")
        )
Fig1E

#Perform pairwise Fisher's exact tests for head and tail offspring separately
library(RVAideMemoire)

#stats for heads
cont.tab_heads<-as.matrix(DeathCount%>% filter(Type=="Head") %>% select(NumDead,NumAlive))
Death_ptable_heads<-fisher.multcomp(cont.tab_heads, p.method = "bonferroni")
Death_ptable_heads

#stats for tails
cont.tab_tails<-as.matrix(DeathCount%>% filter(Type=="Tail") %>% select(NumDead,NumAlive))
Death_ptable_tails<-fisher.multcomp(cont.tab_tails, p.method = "bonferroni")
Death_ptable_tails
```
