I used RSEM to align and quantify RNAseq, then used edgeR to get differentially expressed genes with someting like the following:
library("edgeR");
d0<-read.delim("mydata.txt", row.names="gene_id");
d<-d0[c("c1","c2","t1","t2")];
attach(d);
group<-factor(c(1,1,2,2));
dgeList <- DGEList(counts=d, group=group);
normFac <- calcNormFactors(dgeList); # ERCC normalization
dsgn <- model.matrix(~group);
disp <- estimateDisp(normFac, dsgn);
qlFit <- glmQLFit(disp, dsgn);
qlfTest <- glmQLFTest(qlFit, coef=2);
topTags(qlfTest);
# do GO analysis
go <- goana(qlfTest, specifes="Hs");
The code works up to "topTags(qlfTest)", but stopped at the last line with the following error:
Error in goana.default(qlfTest, specifes = "Hs", de = list(Up = c("A1BG", :
No genes found in universe
Any ideas?
Thanks!
•
written
by
moxu
2 answers
You need to convert your gene identifiers to Entrez gene identifiers. E.g KRT7 is 3855.
I suggest BioMart from Ensembl or https://github.com/stephenturner/annotables
You need to perform a join using the required table (e.g. grch38) and your dataframe of differential expressed genes (qlfTest). The dplyr library is great for jobs like this. A few examples are also included on the github page and http://www.gettinggeneticsdone.com/2015/11/annotables-convert-gene-ids.html.
hello Is there a way to do the goanna and topGo with Ensenmbl or external_gene_name? I lose a lot of genes when I do it with entrez Thanks
I believe your problem lies in that you are using the wrong objects for the goana function.
You cannot just feed it the qflTest and expected it to work. You need to give it a list of genes.
degs <- topTags(qlfTest)[[1]]
go <- goana(degs, species="Hs");
That should work. Ideally in this case you want to feed it:
de = Your differentially expressed genes (up and/or down)
universe = All expressed genes (all row.names of counts/dge object) or alternatively All annotated genes
According to the help information of goana:
de: a vector of Entrez Gene IDs, or a list of such vectors, or an "MArrayLM" fit object.
Yes, you are right.
The problem now is that "topGO(go)" invariantly (I have other differential gene expression analyses) returns the same list of GO pathways as the following:
> go <- goana(as.vector(degs), species="Hs");
> topGO(go)
Term
GO:0000002 mitochondrial genome maintenance
GO:0000003 reproduction
GO:0000009 alpha-1,6-mannosyltransferase activity
GO:0000010 trans-hexaprenyltranstransferase activity
GO:0000012 single strand break repair
GO:0000014 single-stranded DNA endodeoxyribonuclease activity
GO:0000015 phosphopyruvate hydratase complex
GO:0000016 lactase activity
GO:0000018 regulation of DNA recombination
Certainly something is wrong.
Log in to answer this question.
-
Bulk rna sequencing
written by Rosie •
-
Should I subset my data before running EstimateDisp() in RNA-seq analysis ?
written by Basti •
-
EdgeR compare multiple groups in a 1 vs other comparison?
written by shanasabri •
-
Combine multiple groups on comparison by edgeR
written by camelest •
-
Find most upregulated genes in one library?
written by williams.damon26 •
-
How to perform multiple comparisons using edgeR?
written by silas008
-
Too few differentially expressed genes identified by edgeR
written by afli
-
[edgeR usage] Comparing categories of fewer input variables for differential gene expression
written by moxu
-
edgeR Gene Ontology Analysis always returns the same list of pathways
written by moxu
-
Which model to choose in edgeR analysis?
written by biolab
show us a snippet of mydata.txt
Also it should be
species = "Hs"not specifes.Thanks for the correction. However, with "species", I got the same error.
And what's the output of
topTags(qlfTest)?