.Rhistory

library(Hmisc)
?select
?select
?setnames
?rowSums
Matrix::rowSums
?Matrix::rowSums
?select
library(Hmisc)
?select
?get
?names
?aes
?factor
?unique
library(ensembldb)
browseVignettes("ensembldb")
?lapply
?element_text
?duplicted
?duplicated
?length
?sapply
?as.character
?matrix
?unique
conflicts()
conflicts()
library(BiocGenerics)
conflicts()
?geom_bar_repel
library(ggrepel)
?geom_bar_repel
library(ggplot2)
geom_bar_repel
?geom_bar_repel
?get
?select
?scales
?percent
?unit
library(ggplot2)
?unit
class("hello")
test <- list(1,2,3)
test
for(j in list){
for(j in list){
print(j)
}
for(j in test){
print(j)
}
test <- c(1,2,3)
for(j in test){
print(j)
}
test <- c(4,5,2)
for(a in test){}
for(a in test){
print(a)
}
for(a in seq_along(test)){}
for(a in seq_along(test)){
print(a,names(test)[a])
}
for(a in seq_along(test)){
print(a)
}
for(a in seq_along(test)){
print(paste(a,names(test)[a]))
}
for(a in seq_along(test)){
print(paste(a,names(test)[[a]))
for(a in seq_along(test)){
print(paste(a,names(test)[[a]]))
}
for(a in seq_along(test)){
print(paste0(a,names(test)[[a]]))
}
for(a in seq_along(test)){
print(paste0(a,test[[a]]))
}
?seq_along
for(a in seq_along(test)){
print(paste0(a,test[a]))
}
test[2]
c(10,15,20,30)
a <- c(10,15,20,30)
a[a<4]
a[a<30]
library(glue)
test<-"hi"
paste("hello",test)
paste0("hello",test)
paste0("hello",test,"/")
getwd()
library(qs)
data<-qs::qread("../Desktop/allen_Endo.qs")
data
data
colData(data)
colnames(colData(data))
colData(data)["cell_type"]
getwd()
setwd("../Desktop/")
library(qs)
data <- qs::qread("allen_Endo.qs")
data
data
data[, colData(data)[["cluster_celltype"]]=="Micro"]
data[, colData(data)[["cluster_celltype"]]=="asdf"]
file.path()
file.path("test","ok")
getwd()
# IMPORT PACKAGES
library(devtools)
library(roxygen2)
# DIRECTORY FOR PACKAGE
setwd("C:/Users/salma/Power_Analysis_package")
## TESTING ACTUAL PACKAGE
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install()
devtools::check()
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install()
devtools::check()
library(qs)
getwd()
data <- qs::qread("vignettes/data/sce.qs")
data
qsave(data,"sce.qs")
getwd()
qsave(data,"../Desktop/sce.qs")
## TESTING ACTUAL PACKAGE
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install()
devtools::check()
devtools::check()
devtools::build(vignettes=TRUE)
library(poweranalysis)
vignette(poweranalysis)
vignette("poweranalysis")
vignette(package="poweranalysis")
utils::browseVignettes("poweranalysis")
length(list(1,2,3,4))
length(list(1,2,3,4))*2
file.path("test/ow","to/why")
x<-"hi"
file.path(x,"hello","/no/ok")
file.path(x,"why","no")
?file.path
file.path(1,2,3,4)
file.path("output/path/to/file","filename",0.01)
a_0.01 <- 3
a_0.01
5%%3
7%2
7%%2
isOdd <- function(n){}
isOdd <- function(n){
if(n%%2!=0){
return("Odd!")
}else{
return("Not Odd")
}
}
isOdd(434)
isOdd <- function(n){
if(n%%2!=0){
return("Odd!")
}else{
return("Not Odd!")
}
}
isOdd <- function(n){
if(n%%2!=0){
return(paste0(n,"is odd!"))
}else{
return(paste0(n,"is not odd!"))
}
}
isOdd(434)
isOdd <- function(n){
if(n%%2!=0){
return(paste0(n," is odd!"))
}else{
return(paste0(n," is not odd!"))
}
}
isOdd(434)
test <- c(0.1,0.4,-0.1)
for(p in test){}
getwd()
data
colnames(colData(data))
DEout <- DGE_analysis(data,~sex,"sample_id","cluster_celltype")
power_analysis::DEout <- DGE_analysis(data,~sex,"sample_id","cluster_celltype")
DEout <- power_analysis::DGE_analysis(data,~sex,"sample_id","cluster_celltype")
library(power_analysis)
DEout <- poweranalysis::DGE_analysis(data,~sex,"sample_id","cluster_celltype")
DEout
DEout
sex_chromosome_DEGs <- function(all_genes,
ensemblID_colname="name"){
# check input parameter is fine
if(ensemblID_colname!="name"){
if(!is.character(ensemblID_colname)){
stop("Error: ensemblID_colname should be a string specifying the column name for the column containing ensembl IDs.")
}
}
if(class(all_genes)!="list"){
stop("Error: all_genes should be a list")
}
column_exists <- sapply(all_genes, function(x) ensemblID_colname %in% names(x))
if(!all(column_exists)){
stop("Error: every DGE analysis output in all_genes should contain the column with the name specified by ensemblID_colname")
}
# combine data for all celltypes into one dataframe
combn_genes <- rbindlist(all_genes,idcol="celltype")
# add symbol
mart <- useDataset("hsapiens_gene_ensembl", useMart("ensembl"))
genes <- unique(as.character(combn_genes$name))
gene_IDs <- getBM(filters = "ensembl_gene_id",
attributes = c("ensembl_gene_id","hgnc_symbol","chromosome_name"),
values = genes, mart = mart, useCache = FALSE)
gene_IDs <- as.data.table(gene_IDs)
setnames(gene_IDs,"ensembl_gene_id",ensemblID_colname)
# remove any duplicates from reference set - two names for one ENSEMBL ID
gene_IDs <- unique(gene_IDs,by=ensemblID_colname)
setkey(combn_genes,name)
# append gene names
combn_genes[, gene_name := gene_IDs[combn_genes, on=.(name), x.hgnc_symbol]]
# append chromosome
combn_genes[, chromosome := gene_IDs[combn_genes, on=.(name), x.chromosome_name]]
setorder(combn_genes,celltype,adj_pval,logFC)
# include only genes on the sex chromosomes
sexgenes <- combn_genes[combn_genes$chromosome == "X" | combn_genes$chromosome == "Y"]
# remove unnecessary columns and reshape this as appropriate
sexgenes <- split(sexgenes,sexgenes$celltype)
sexgenes <- sexgenes[[names(sexgenes)]][, c("logFC","logCPM","LR","PValue","adj_pvalue","name")]
# output DGE analysis subset to sex chromosome genes
return(sexgenes)
}
DEout$celltype_DEGs
DEout$celltype_DEGs$Micro
DEout$celltype_DEGs
sex_chromosome_DEGs(DEout$celltype_DEGs$Micro)
sex_chromosome_DEGs(DEout$celltype_DEGs$Micro$name)
sex_chromosome_DEGs(list(DEout$celltype_DEGs$Micro$name))
sex_chromosome_DEGs(DEout$celltype_all_genes)
library(data.table)
sex_chromosome_DEGs(DEout$celltype_all_genes)
library(biomaRt)
sex_chromosome_DEGs(DEout$celltype_all_genes)
sex_chromosome_DEGs(DEout$celltype_all_genes)
data
colnames(colData(data))
colData(data)["cluster_celltype"]
data[["Micro"]]
SCEs <- list(data)
name(data)
names(data)
names(SCEs)
for(study in names(SCEs)){
print(study)
}
for(study in SCEs){
print(study)
}
celltype_correspondence <- list(
"Microglia" = list("micro", "Micro", "microglia"),
"Astrocytes" = list("astro", "Astro", "astrocytes")
)
celltype_correspondence[["Microglia"]]
celltype_correspondence[["Microglia"]][[1]]
celltype_correspondence[["Microglia"]][[2]]
data[, colData(data)[["cluster_celltype"]] == "Micro"]
data[, colData(data)[["cluster_celltype"]] == "poop"]
class(data)
class(data)[[1]]
DEout
names(DEouts)
names(DEou)
names(DEout)
list(1,1,2)
data
df
DEout
getwd()
setwd(~/Desktop)
setwd("../Desktop/")
save(DEout,file="DEout.RData")
if(3==4)
print("no")
else
print("ok")
if(3==4)
print("why")
if(4==2)
if(7==9)
print("no")
if(4==2)
if(7==9)
print("no")
x=TRUE
if(isTrue(x))
print("ok")
if(isTRUE(x))
print("ok")
x=F
if(isTRUE(x))
print("ok")
x=T
if(isTRUE(x))
print("ok")
if(x){}
if(x){
print("why")
}
x
x=T
x
DEout
DEout2<-DEout
test <- list(DEout,DEout2)
names(test)
test
for(i in names(test)){
print(i)
}
test
test <- list(DEout,DEout2)
de1 <- DEout$celltype_all_genes
de2 <- DEout2$celltype_all_genes
test <- list(de1,de2)
for(i in names(test)){
print(i)
}
for(i in names(test)){
print(i)
}
names(de1)
names(de2)
test[[1]]
test <- list(DEout,DEout2)
names(test[[1]]$celltype_all_genes)
x <- list(1,2,3)
is.vector(x)
test <- list(DEout,DEout2)
for(i in seq_along(DEouts)){
# get corresponding cell type names for current study
celltype_name <- celltypes_list[[i]]
# redefine DEouts so each element only contains study/celltype
allstudies[[i]] <- DEouts[[i]]$celltype_all_genes[[celltype_name]]
}
# variable to redefine DEouts
allstudies <- list()
# for each study, select data corresponding only to celltype
for(i in seq_along(test)){
# get corresponding cell type names for current study
#celltype_name <- celltypes_list[[i]]
# redefine DEouts so each element only contains study/celltype
allstudies[[i]] <- test[[i]]$celltype_all_genes[["Micro"]]
}
allstudies
names(allstudies) <- c("data1","data2")
allstudies
# reshape data so "dataset" is now a variable
allstudies_dt <- rbindlist(allstudies,idcol="dataset")
allstudies_dt
unique(allstudies_dt$dataset)
allGenes <- list()
for(j in 1:length(data_names)){
allGenes[[j]] <- allstudies_dt[dataset==data_names[[j]],name]
}
data_names <- c("data1","data2")
allGenes <- list()
for(j in 1:length(data_names)){
allGenes[[j]] <- allstudies_dt[dataset==data_names[[j]],name]
}
allGenes
allGenes[[1]]
allGenes[[2]]
allGenes[[3]]
# get set of genes in all studies
shared_genes <- Reduce(intersect,allGenes)
shared_genes
setkey(allstudies_dt,name)
allstudies_dt <- allstudies_dt[shared_genes]
# select DEGs by specifying only significant genes
genes <- allstudies_dt[dataset==dataset_name & PValue<0.05, name]
# filter to just these genes
allstudies_dt <- allstudies_dt[name %in% genes,]
dataset_name = "data1"
# select DEGs by specifying only significant genes
genes <- allstudies_dt[dataset==dataset_name & PValue<0.05, name]
# filter to just these genes
allstudies_dt <- allstudies_dt[name %in% genes,]
# make matrix for corr() - cols will be [datset, name, logFC]
mat_lfc <- allstudies_dt[,.(dataset,name,logFC)]
# reshape so cols now are [(gene) name, logFC.dataset1, logFC.dataset2,...,logFC.datasetN] (N being length(names(DEouts)))
mat_lfc <-
reshape(mat_lfc, idvar = "name", timevar = "dataset",
direction = "wide")
# remove logFC. from name (now cols are just [name, dataset1, dataset2,...,datasetN] with logFC values in each column)
colnames(mat_lfc)<-
c(colnames(mat_lfc)[1],
substr(colnames(mat_lfc[,2:ncol(mat_lfc)]),
7,nchar(colnames(mat_lfc[,2:ncol(mat_lfc)]))))
# remove NA rows // convert NA to 0s
mat_lfc <- mat_lfc[complete.cases(mat_lfc), ]
# get number of DEGs for this celltype
num_genes <- dim(mat_lfc)[[1]]
print(paste0("Number of ",celltype," genes at given p value is ",num_genes))
genes <- mat_lfc$name # (only updates if we remove NA genes, else is the same as "genes" above)
# get correlation matrix
corr_lfc <- cor(mat_lfc[,2:ncol(mat_lfc)],method = "spearman") #### could also try "pearson"...
# get number of DEGs for this celltype
num_genes <- dim(mat_lfc)[[1]]
print(paste0("Number of ",celltype," genes at given p value is ",num_genes))
genes <- mat_lfc$name # (only updates if we remove NA genes, else is the same as "genes" above)
celltype="Micro"
# get number of DEGs for this celltype
num_genes <- dim(mat_lfc)[[1]]
print(paste0("Number of ",celltype," genes at given p value is ",num_genes))
genes <- mat_lfc$name # (only updates if we remove NA genes, else is the same as "genes" above)
# get correlation matrix
corr_lfc <- cor(mat_lfc[,2:ncol(mat_lfc)],method = "spearman") #### could also try "pearson"...
corr_lfc
celltype_correspondence <- list(
"Microglia" = list("micro", "Micro"),
"Astrocytes" = list("astro", "Astro")
)
celltype_correspondence [["Microglia"]]
celltype_correspondence [["Microglia"]][[1]]
for(celltype in celltype_correspondence){
print(celltype)
}
celltype_correspondence
names(celltype_correspondence)
celltype_correspondence[["Microglia"]]
setwd("../Power_Analysis_package/")
## TESTING ACTUAL PACKAGE
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install(build_vignettes=TRUE)
## TESTING ACTUAL PACKAGE
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install(build_vignettes=TRUE)
## TESTING ACTUAL PACKAGE
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install(build_vignettes=TRUE)
## TESTING ACTUAL PACKAGE
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install(build_vignettes=TRUE)
devtools::check()
devtools::load_all()
devtools::document()
devtools::check()
sce<-qs::qread("vignettes/data/sce.qs")
qsave(sce,file=sce.qs)
qsave(sce,file="sce.qs")
qsave(sce,file="vignettes/data/sce.qs")
## TESTING ACTUAL PACKAGE
devtools::load_all()
devtools::document()
devtools::check()
## TESTING ACTUAL PACKAGE
devtools::load_all()
devtools::document()
devtools::check()
## TESTING ACTUAL PACKAGE
devtools::document()
devtools::load_all()
devtools::check()
getwd()
Sys.setFileTime(".", Sys.time())
.rs.restartR()
getwd()
getwd()
## TESTING ACTUAL PACKAGE
devtools::document()
devtools::load_all()
devtools::check()
## TESTING ACTUAL PACKAGE
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install()
devtools::document()
#usethis::use_vignette("poweranalysis")
devtools::load_all()
devtools::install(use_vignettes = TRUE)
devtools::install(use_vignette = TRUE)
devtools::check()
system("R CMD build .")
system("R CMD check --as-cran mypackage_*.tar.gz")