goseq with non native organism

Question

[alyyha]

I´m using goseq package after DESeq2 analysis and I am totally new in R

When I would like to test the pwf of my DE genes, I get an error:

if (!requireNamespace("BiocManager", quietly = TRUE))
install.packages("BiocManager")

BiocManager::install("TxDb.Athaliana.BioMart.plantsmart28")

library(TxDb.Athaliana.BioMart.plantsmart28)
txdb<-TxDb.Athaliana.BioMart.plantsmart28
txsByGene=transcriptsBy(txdb,"gene")
lengthData=median(width(txsByGene))



assayed.genes <- rownames(res)
de.genes <- rownames(res)[ which(res$padj < 0.05)]
genes<-as.integer(assayed.genes%in%de.genes)
names(genes)<-row.names(res$padj)

genes
    0     1 
29357  3476 

pwf<-nullp(genes,txdb,assayed.genes,lengthData)

Error in nullp(genes, txdb, assayed.genes, lengthData) : 
bias.data vector must have the same length as DEgenes vector!

length(lengthData)=33602

length(genes)=32833

DEgenes vector in the script is genes and contains all Arabidopsis genes (32833). I did not filter for counts.

bias.data vector in the script is lengthData and contains median length values of all transcripts (33602).

Does someone knows how to solve this error? Thank you very much in advance!

Answer 1

[ATpoint]

You have 32833 genes and 33602 elements in the bias vector, see the problem? You must have one value in the bias bector and it must have the exact same order as the gene list.

Code suggestion:

assayed.genes <- rownames(res)
genes         <- as.integer(res$padj < 0.05)
biasdata      <- lengthData[match(assayed.genes, names(txsByGene))]
nP            <- nullp(DEgenes = genes, bias.data = biasdata)

[alyyha]

Yes I see this problem, thank you. I am just wondering how to fix it because I don´t know how to do it.

[ATpoint]

Just added a suggestion.

[alyyha]

Thank you very much! I'll try

[ATpoint]

Also added a command for nullp. The DEgenes is a vector of same length as biasdata and it has a 1 if the gene is DE and a 0 if not.

[alyyha]

It worked. Danke schon!!

[ATpoint]

Glad to help! If the answer solved the issue please accept it to provide closure (using the checkmark button).