06_annotationBlocksCount_CGI.R 1.28 KB
Newer Older
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
library(rtracklayer)
library(stringr)

# CpG island file, bed to GRanges conversion -----------
cgi <- import("ucsc_Sscrofa11.2_cpgislands.bed.gz", format = "BED")
seqlevels(cgi) <- str_remove(seqlevels(cgi), "chr")
saveRDS(cgi, file = "ucsc_Sscrofa11.2_cpgislands.RDS")

# sarray script -----------

library(dplyr)
library(purrr)
library(stringr)

prefix <- "/work2/genphyse/genepi/data/porc/medpseq/ROSEpigs_processed2/results/bwa/mergedLibrary/"
prefix <- "/home/gdevailly/mnt/genotoul_genepi/data/porc/medpseq/ROSEpigs_processed2/results/bwa/mergedLibrary/"

bam_data <- tibble(
    path = list.files(prefix, pattern = "*.bam$"),
    id = map_chr(strsplit(path, ".", fixed = TRUE), 1),
    type = str_sub(id, end = -4)
)

commands <- paste0(
    "./annotationBlocksCount.R -i ",
    prefix,
    bam_data$path,
    " -w /work2/genphyse/genepi/guillaume/rosepigs/ucsc_Sscrofa11.2_cpgislands.RDS -o rds_cgi/",
    bam_data$id, "_cgi.RDS"
)

write.table(
    commands,
    file = "annotationBlocksCount_CGI_sarray.sh",
    quote = F, row.names = F, col.names = F
)

# chmod +x annotationBlocksCount_CGI_sarray.sh
# module load compiler/gcc-7.2.0
# module load module load system/R-4.0.2_gcc-7.2.0
# sarray -J abc -o %j.out -e %j.err -t 01:00:00 --mem=32G --mail-type=FAIL annotationBlocksCount_CGI_sarray.sh