-
Notifications
You must be signed in to change notification settings - Fork 2
/
scrap.rmd
69 lines (40 loc) · 2.14 KB
/
scrap.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
# Scrap
Misc code removed fromt he main manuscript that may still be of use later.
### Taxa
The primary intent of comparing taxa across matrices was to validate our taxon name reconciliation across studies.
We first considered pairwise similarity between the same species from different matrices in different studies.
```{r taxon_mapping_across_studies, eval=FALSE}
Dtaxaxman =
read_tsv("../reconciliation/blast/graphs/cross_mscript_taxa.tsv") %>%
arrange( mscript1, dataset1, taxon1, mscript2, dataset2, taxon2 )
Dtaxaxman_raw_count = nrow(Dtaxaxman)
Dtaxaxman %<>%
distinct( mscript1, dataset1, taxon1, mscript2, dataset2, taxon2, .keep_all = TRUE )
```
`r Dtaxaxman_raw_count` raw pairwise comparisons passed the similarity thresholds under these sampling criteria. These included multiple hits for the same sequence pairs from different regions. When only a single hit is retained, there are `r nrow(Dtaxaxman)` retained sequence pairs.
We next considered taxon mapping more generally.
```{r taxon_mapping, eval=FALSE}
Dtaxa =
read_tsv("../blast/graphs/taxon_graph.tsv") %>%
arrange( mscript1, dataset1, taxon1, mscript2, dataset2, taxon2 )
# Set aside species correspondense to same species in same matrix in same matrix
Dtaxa_self =
Dtaxa %>%
filter( mscript1 == mscript2 & dataset1==dataset2 & taxon1==taxon2 )
# Remove species correspondense to same species in same matrix in same study
Dtaxa %<>%
filter( !(mscript1 == mscript2 & dataset1==dataset2 & taxon1==taxon2) )
n_not_self = nrow(Dtaxa)
# Remove hits in same matrix in same study
Dtaxa %<>%
filter( !(mscript1==mscript2 & dataset1==dataset2))
# Remove hits same study, presume that taxon resolution is OK within study.
Dtaxa %<>%
filter( !(mscript1==mscript2))
ggplot(Dtaxa) + geom_density(aes(count)) + xlim(0,50)
ggplot(Dtaxa_self) + geom_density(aes(count))
# Top of this list are candidates for renaming because they are very similar but have different names
Dtaxa %>% filter(taxon1!=taxon2) %>% arrange( desc(count))
# Top of this list are candidates for renaming because they are very different but have same names
Dtaxa %>% filter(taxon1==taxon2) %>% arrange( count )
```