/
paristools.Rmd
executable file
·82 lines (56 loc) · 1.4 KB
/
paristools.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
---
title: "Introduction to paristools"
author: "Zhuoer Dong"
date: "`r Sys.Date()`"
output: prettydoc::html_pretty
vignette: >
%\VignetteIndexEntry{Introduction to paristools}
%\VignetteEngine{knitr::rmarkdown}
%\VignetteEncoding{UTF-8}
---
```{r setup, include = FALSE}
knitr::opts_chunk$set(collapse = TRUE, comment = "#>")
```
```{r message=FALSE}
library(paristools)
library(ggplot2)
```
# read files
```{r}
sam_file <- system.file('extdata', 'Neat1_1.Aligend_trunc.sam', package = 'paristools')
mapping <- read_sam(sam_file)
mapping
```
```{r}
duplexgroup_file <- system.file('extdata', 'Neat1_1.duplexgroup', package = 'paristools');
duplex_group <- read_duplexgroup(duplexgroup_file)
duplex_group
```
refer to `read_bed()`
# calculate coverage
## reference coverage
```{r}
loc_df <- mapping |> sam_to_loc_df()
loc_df
ref_coverage <- loc_df |> cal_coverage()
ref_coverage
```
```{r}
ggplot(ref_coverage) +
geom_col(aes(pos, n_reads)) +
facet_wrap(~strand, ncol = 1) +
labs(x = 'position', y = 'number of reads')
```
## duplex group segment average coverage
```{r}
dg_coverage <- get_dg_coverage(duplex_group, ref_coverage)
dg_coverage
```
```{r}
dg_coverage |> dplyr::select(-strand) |> tidyr::spread('pair', 'n_reads') |>
ggplot() + geom_point(aes(left, right))
```
# others
```{r}
mapping |> dplyr::mutate(strand = purrr::map_chr(FLAG, paristools:::get_strand))
```