-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy path03-G3.rmd
83 lines (60 loc) · 2.98 KB
/
03-G3.rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
# IGHV1-24 - G3
```{r}
source("functions.R")
```
```{r, echo=FALSE}
load("data.rda")
chain = "IGH"
func <- data.frame(allele = names(vgerms[[chain]]), functionality = !grepl("(ORF|P)",sapply(seqinr::getAnnot(vgerms_full[[chain]]), function(x) unlist(strsplit(x,"[|]"))[4])), sign = sapply(seqinr::getAnnot(vgerms_full[[chain]]), function(x) unlist(strsplit(x,"[|]"))[4]), stringsAsFactors = F)
mat <- mat_list$IGH$functional$nonsingle$all$`318`
load("data_frac_new2.rda")
data <- setDT(data_frac$IGH$functional$nonsingle$all$`318`$complete$`95`)
data[,v_call:=paste0(v_gene,"*",v_allele)]
load("alleles_dbs.rda")
allele_db <- alleles_dbs$IGH$functional$nonsingle$all$`318`$complete$`95`
allele_db <- allele_db %>% dplyr::rowwise() %>% dplyr::mutate(gene = alakazam::getGene(or_allele, strip_d = F, omit_nl = F), group = strsplit(gsub(gene, "", new_allele),"[*]")[[1]][1], gene_group = alakazam::getGene(new_allele, strip_d = F, omit_nl = F))
load("functional_groups.rda")
func_groups <- functional_groups$IGH$functional$nonsingle$all$`318`$complete$`95`
cols <- c("#FAAB18", "#1380A1","#990000", "#588300")
pal <- cols %>%
newpal(names = c("orangy", "bluish", "redish", "greeny"))
edit_links <- readLines("edit_links.txt")
share_links <- readLines("share_links.txt")
```
```{r,echo=FALSE}
g_group = 'IGHV1-24G3'
group = names(func_groups)[func_groups==g_group]
gr <- allele_db %>% filter(gene_group == g_group) %>% pull(group) %>% unique()
g <- allele_db %>% filter(gene_group == g_group) %>% pull(gene) %>% unique()
```
## Allele appearnce
The group of `r group` includes `r length(grep(g_group,allele_db$new_allele,value=T))` alleles, `r sum(func$functionality[func$allele %in% allele_db$or_allele[grep(g_group,allele_db$new_allele)]])` out of the alleles are functional.
For each allele we counted the number of appearances across the population, any appearance was considered valid.
```{r}
allele_appearance(data, g_group, allele_db)
```
## Sequence depth
To examine the potential cutoff we observed the sequence depth for each allele
```{r}
tagList(sequence_depth(data, g_group, allele_db))
```
## Observations
For this group we've seen a single allele, hence the general summary statistics do not apply here and were omitted.
This section is editable by clicking on the edit button below. To refresh the section click on the refresh button
You can access the file also from [here](`r edit_links[grep(g_group,edit_links)]`){target="_blank"}
```{r, message=F, echo=F, warning=F}
knitr::include_url(paste0("https://peresay.shinyapps.io/G2_group/?group=",g_group))
```
## Conclusions
Because we observed a single allele in this group we set the cutoff to $0.01%$ absolute threshold (the relative repertoire frequency). The plot below shows the frequency of the group above the threshold line.
```{r, out.height= "100%", out.width="100%"}
m <- list(
l = 50,
r = 50,
b = 100,
t = 100,
pad = 0.5
)
heatmap_alleles(data, g_group, allele_db, func)%>%
layout(autosize = F, width = 800, height = 1300)
```