Hi,
I am a novice for R and bioinfomatics. I am trying to run DESeq2 with my raw count table (.txt). I have been trying to follow the beginner's guide for the DESeq2 package, but it is still hard to understand because my experimental condition is different from the example. Would you please help me to figure it out how to make a DESeqDataset and run? Thank you in advance for your time.
Here is my experiemtal design.
- Treatment: Control vs Treated
- Animal group: 4 genotypes in duplicate (genotype1-1, genotype1-2, genotype2-1, genotype2-2, genotype3-1, genotype3-2, genotype4-1, genotype4-2)
- Model: Treatment, Animal group, Aniaml * Treatment (interaction of animal and treatment).
I want to run DESeq2 with 2-way ANOVA-like frame work (2 treatment X 4 genotypes).
Jin
Which animal groups are treated and untreated?
Each treatment has 4 genotypes.
In other words, Control (genotype1-1, genotype1-2, genotype2-1, genotype2-2, genotype3-1, genotype3-2, genotype4-1, genotype4-2).
Treated (genotype1-1, genotype1-2, genotype2-1, genotype2-2, genotype3-1, genotype3-2, genotype4-1, genotype4-2).
So if I understand correctly, you have total 16 samples. 8 treated and 8 control, each with 4 genotypes in duplicates.
But again if you post this on bioconductor support forum you will get proper suggestion from developer himself.
Hi poisonAlien,
Thank you very much for your help. I will try this and post on the forum as your suggestion as well.
Jin
I tried it, and got a following Error like below.
> dds = DESeqDataSetFromMatrix(countData = countMatrix, colData = coldata, design = ~ group + treatment + group:treatment)
Error in as.matrix(countData) :
error in evaluating the argument 'x' in selecting a method for function 'as.matrix': Error: object 'countMatrix' not found
Calls: DESeqDataSetFromMatrix -> as.matrix
countMatrix is your counts table. (as you mentioned in your question - raw count table (.txt).) You should read it into R before doing this.
Assuming you have named your table as 'count_table.txt' and first column is your gene names (or IDs), followed by rest of the 16 columns, organised according to sample names in above 'coldata'
countMatrix = read.delim("count_table.txt",sep='\t',row.names = 1)
then try running above code.
There was no error, but give some message below,
> dds = DESeqDataSetFromMatrix(countData = countMatrix, colData = coldata, design = ~ group + treatment + group:treatment)
it appears that the last variable in the design formula, 'treatment', has a factor level, 'control', which is not the base level. we recommend you use factor(...,levels=...) or relevel() to set this as the base level before proceeding. for more information, please see the 'Note on factor levels' in vignette('DESeq2').
Ahh ! I The warning because reference level is treated here. See the edit above.
It works! Thank you so much for your help!
I have got the result below,
Does this mean that the genotype4_treated is significantly different from others?
log2 fold change (MAP): groupgt4.treatmenttreated
Wald test p-value: groupgt4.treatmenttreated
DataFrame with 39595 rows and 6 columns
baseMean log2FoldChange lfcSE stat pvalue padj
<numeric> <numeric> <numeric> <numeric> <numeric> <numeric>
unigene22398972 189645.38 -0.1773434 0.2744156 -0.6462584 0.518111997 0.999899
unigene22407675 179211.24 0.1429904 0.1045016 1.3683091 0.171215344 0.999899
unigene22402019 72765.83 -0.1883720 0.2830502 -0.6655075 0.505725927 0.999899
unigene22387256 70438.23 0.3585405 0.1246933 2.8753784 0.004035435 0.999899
unigene22410725 29215.47 0.2276120 0.1606212 1.4170736 0.156461409 0.999899
... ... ... ... ... ... ...
unigene22407834 0.07103343 -5.161589e-05 0.01173293 -0.004399233 0.9964899 NA
unigene22405341 0.07103343 -5.161589e-05 0.01173293 -0.004399233 0.9964899 NA
unigene22407238 0.07103343 -5.161589e-05 0.01173293 -0.004399233 0.9964899 NA
unigene22407247 0.07103343 -5.161589e-05 0.01173293 -0.004399233 0.9964899 NA
unigene22397119 0.07103343 -5.161589e-05 0.01173293 -0.004399233 0.9964899 NA
I can't say much about your results. But regarding the above error, your first column must start with counts. Here your first column is ID's, which should be rowanames. While reading the table, use
row.names = 1
, which makes first columns as rowanmes. Anyways here is it.Thank you so much for your help again. I appreciate it.
Hi, I am just bit confused.. I run the following script on R version 3.2 and everything was fine:
however, when I updated R into 3.3 I start getting the following error message:
I tried to follow many solutions but non of them worked .
can anyone help
thanks much