Entering edit mode
19 months ago
yuhangzhu199512
•
0
Hi, I am using clusterProfiler to do an enrichment analysis. However, I get the warning ' No gene can be mapped.... --> Expected input gene ID: ND6,ND4,COX3,ND1,ND5,ATP6 --> return NULL...'. I am studying a non-model organism, the 'Vulpes lagopus'. My script is as follows,
mycounts<-read.csv("-B5_Counts-Arctic-fat.csv",row.names = 1)
head(mycounts)
dim(mycounts)
library(dplyr)
library(DESeq2)
mycounts <- mycounts %>%mutate_if(is.factor, ~as.numeric(as.character(.)))
mymeta<-read.csv("-B5_mymeta.csv",stringsAsFactors = T)
mycountsGroups <- mymeta$Group
colData <- data.frame(row.names=colnames(mycounts), group_list=mycountsGroups)
dds <- DESeqDataSetFromMatrix(countData = mycounts, colData = colData, design = ~group_list, tidy = F)
keep <- rowSums(counts(dds)>10) >4
dds.keep <- dds[keep, ]
dim(dds.keep)
degobj <- DESeq(dds.keep);degobj
countdds <- counts(degobj, normalized = T)
deg.res1 <- results(degobj, contrast=c("group_list","HFB","LFB"), lfcThreshold = 0, pAdjustMethod = 'fdr')
dim(deg.res1)
apply(deg.res1,2, function(x) sum(is.na(x)))
deg.res1$padj[is.na(deg.res1$padj)] <- 1
table(is.na(deg.res1$padj))
summary(deg.res1)
library(dplyr)
res_1<-data.frame(deg.res1)
res_1 %>% mutate(group = case_when(log2FoldChange >0.5 & padj <= 0.05 ~ "UP",log2FoldChange <-0.5 & padj <= 0.05 ~ "DOWN",TRUE ~ "NOT_CHANGE")) -> res_2
table(res_2$group)
upsig <- res_2[grep(pattern = "UP",res_2[,7]),]
downsig <- res_2[grep(pattern = "DOWN",res_2[,7]),]
deg <- rbind(upsig,downsig)
BiocManager::install("AnnotationHub")
BiocManager::install("biomaRt")
# load package
library(AnnotationHub)
library(biomaRt)
# make a orgDb
hub <- AnnotationHub::AnnotationHub()
query(hub, "Vulpes lagopus")
Vulpes.OrgDb <- hub[["AH96456"]]
columns(Vulpes.OrgDb)
library(clusterProfiler)
df2 <- bitr(row.names(upsig), fromType = "SYMBOL",
toType = c("ENTREZID"),
OrgDb = Vulpes.OrgDb)
genelist <- df2$ENTREZID
enrich.go.BP = enrichGO(gene =genelist,
OrgDb = Vulpes.OrgDb,
keyType = 'ENTREZID',ont= "BP",
pvalueCutoff = 0.05,
pAdjustMethod = "BH",
qvalueCutoff = 0.05,
)
--> No gene can be mapped....
--> Expected input gene ID: 23629751,23629755,23629750,23629748,23629754,23629750
--> return NULL...
cross-posted on bioconductor: link