Deseq2 : Script out of bounds
1
0
Entering edit mode
3.3 years ago
aranyak111 • 0

I am using a code based on Deseq2. One of my goals is to plot a heatmap of data.

heatmap.data <- counts(dds)[topGenes,]

The error I am getting is

Error in counts(dds)[topGenes, ]: subscript out of bounds

the first few line sof my counts(dds) function looks like this.

                      99h1   99h2   99h3   99h4   wth1   wth2
ENSDARG00000000002    243    196    187    117     91     96
ENSDARG00000000018     42     55     53     32     48     48
ENSDARG00000000019     91     91    108     64     95     94
ENSDARG00000000068      3     10     10     10     30     21
ENSDARG00000000069     55     47     43     53     51     30
ENSDARG00000000086     46     26     36     18     37     29
ENSDARG00000000103    301    289    289    199    347    386
ENSDARG00000000151     18     19     17     14     22     19
ENSDARG00000000161     16     17      9     19     10     20
ENSDARG00000000175     10      9     10      6     16     12
ENSDARG00000000183     12      8     15     11      8      9
ENSDARG00000000189     16     17     13     10     13     21
ENSDARG00000000212    227    208    259    234     78     69
ENSDARG00000000229     68     72     95     44     71     64
ENSDARG00000000241     71     92     67     76     88     74
ENSDARG00000000324     11      9      6      2      8      9
ENSDARG00000000370     12      5      7      8      0      5
ENSDARG00000000394    390    356    339    283    313    286
ENSDARG00000000423      0      0      2      2      7      1
ENSDARG00000000442      1      1      0      0      1      1
ENSDARG00000000472     16      8      3      5      7      8
ENSDARG00000000476      2      1      2      4      6      3
ENSDARG00000000489    221    203    169    144     84    114
ENSDARG00000000503    133    118    139     89     91    112
ENSDARG00000000529     31     25     17     26     15     24
ENSDARG00000000540     25     17     17     10     28     19
ENSDARG00000000542     15      9      9      6     15     12

How do I ensure all the elements of the top genes are present in it?

When I try to see 20 top genes in the dataset. it looks like a list of genes

 [1] "6339"  "12416" "1241"  "3025"  "12791" "846"   "15090"
 [8] "6529"  "14564" "4863"  "12777" "1122"  "7454"  "13716"
[15] "5790"  "3328"  "1231"  "13734" "2797"  "9072"

with the column head V1.

I have used both

topGenes <- read.table("E://mir99h50 Cheng data//topGenesresordered.txt",header = TRUE)

and

topGenes <- read.table("E://mir99h50 Cheng data//topGenesresordered.txt",header = FALSE)

to see if the out of bounds error is removed. However it was of no use. I guess the V1 head is causing the issue.

The top genes function has been generated using the above code snippet.

resordered <- res[order(res$padj),]

#Reorder gene list by increasing pAdj
resordered <- as.data.frame(res[order(res$padj),])

#Filter for genes that are differentially expressed with an FDR < 0.01
ii <- which(res$padj < 0.01)
length(ii)

# Use the rownames() function to get the top 20 differentially expressed genes from our results table

topGenes <- rownames(resordered[1:20,])
topGenes

# Get the counts from the DESeqDataSet using the counts() function
heatmap.data <- counts(dds)[topGenes,]
Genomics deseq2 R • 2.3k views
ADD COMMENT
0
Entering edit mode

The error message is pretty clear - you're trying to use an invalid index. Check if all elements of topGenes are present in counts(dds).

ADD REPLY
0
Entering edit mode

I moved your additional content from the answer field to this post. Please in the future use edit to add content, not the answer field, otherwise the thread quickly becomes unorganized.

As for the question, the important detail that is missing is what the content of topGenes is and how it was created. Please comment on existing answers and comments via ADD COMMENT/REPLY.

ADD REPLY
0
Entering edit mode

The information you've added in has only caused more confusion, not explained anything. Is your topGenes obtained from the rownames of resorderedor read from a file with read.table?

What is the command you use when you get the list of numbers you describe as "When I try to see 20 top genes..."?

ADD REPLY

Login before adding your answer.

Traffic: 2340 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6