In byandell/DOread: Read DO data

knitr::opts_chunk$set(echo = TRUE)

suppressPackageStartupMessages({
  library(DOread)
  library(feather)
  library(readxl)
  library(dplyr)
})

datapath <- "~/Documents/Research/attie_alan/DO/"
dirpath <- file.path(datapath, "data", "DerivedData")

if(!dir.exists(file.path(dirpath, "RNAseq")))
  dir.create(file.path(dirpath, "RNAseq"))

# Old data
#load(file.path(dirpath, "RNAseq", "DO381_islet.RData"))
annot.samples$Mouse.ID <- stringr::str_replace(annot.samples$Mouse.ID, "DO-0*", "")
head(annot.samples)
#There are mouse duplicates. Throw them all out.
(dupes <- sort(unique(c(which(duplicated(annot.samples$Mouse.ID)),
                  which(duplicated(annot.samples$Mouse.ID, fromLast = TRUE))))))
if(length(dupes)) {
  annot.samples <- annot.samples[-dupes, ]
  expr.mrna <- expr.mrna[-dupes, ]
}
saveRDS(annot.samples,
        file = file.path(dirpath, "RNAseq", "annot.samples.rds"))

load(file.path(datapath, "AttieDOv2", "DerivedData", "Attie_islet_secr_data_v1.Rdata"))

ls()

head(annot.mrna)

Add covar to annot.mrna.

annot.mrna$covar <- "sex,DOwave"

Save annotations in RDS files.

saveRDS(annot.mrna, file = file.path(dirpath, "RNAseq", "annot.mrna.rds"))

Save rankz.mrna in feather file. First create mouse ID column

rankz.mrna <- data.frame(Mouse.ID = rownames(rankz.mrna), rankz.mrna,
                        stringsAsFactors = FALSE)

feather::write_feather(rankz.mrna, file.path(dirpath, "RNAseq", "expr.mrna.feather"))

Expression Peaks

# Old file
system.time(peaks <- read_excel(file.path(datapath, "rna_seq", "all_suggestive_qtls_mrna.xlsx")))
table(table(peaks$gene_id))

peaks.mrna <- readRDS(file.path(datapath, "AttieDOv2", "DerivedData", "rnaseq_rankz_peaks.rds"))

pheno_tissue <- DOread:::pheno_tissue

Positions in bp, not Mbp.

m <- match(peaks.mrna$lodcolumn, annot.mrna$id)
peaks.mrna$symbol <- annot.mrna$symbol[m]
peaks.mrna$gene_chr <- annot.mrna$chr[m]
peaks.mrna$gene_start <- annot.mrna$start[m]
peaks.mrna$gene_end <- annot.mrna$end[m]
peaks.mrna$pos <- peaks.mrna$pos

tissue <- "Islet"
chrs <- names(map)
peaks.mrna <- peaks.mrna %>%
  mutate(chr = factor(chr, chrs),
         pheno = pheno_tissue(tissue, symbol, lodcolumn),
         longname = pheno,
         output = pheno,
         pheno_group = "Islet.mRNA",
         pheno_type = paste("Islet.mRNA", gene_chr, sep = ".")) %>%
  rename(gene_id = lodcolumn,
         qtl_chr = chr,
         qtl_pos = pos) %>%
  select(-lodindex)

feather::write_feather(peaks.mrna, file.path(dirpath, "RNAseq", "peaks.mrna.feather"))

byandell/DOread documentation built on May 13, 2019, 9:26 a.m.