I am performing differential expression analysis with DESeq2. After running dds<- DESeq(dds)
i wanted to retrive all the up and down regulated genes. i gave the conditionpvalue < 0.05 and log2FoldChnage > 1.5 for up regulated genes. But in the result file also include the genes which have pvalue > 0.05. iam pasting my example dataset here with the result. I wnted to know that which command will be helpfull here for retriving genes which satisfy the above mentioned two condition.
Thank you in advance.
Code should be placed in three backticks as shown below
> dds <- makeExampleDESeqDataSet()
> dds
class: DESeqDataSet
dim: 1000 12
metadata(1): version
assays(1): counts
rownames(1000): gene1 gene2 ... gene999 gene1000
rowData names(3): trueIntercept trueBeta trueDisp
colnames(12): sample1 sample2 ... sample11 sample12
colData names(1): condition
> dds <- DESeq(dds)
estimating size factors
estimating dispersions
gene-wise dispersion estimates
mean-dispersion relationship
final dispersion estimates
fitting model and testing
> res <- results(dds, alpha =0.05)
> res
log2 fold change (MLE): condition B vs A
Wald test p-value: condition B vs A
DataFrame with 1000 rows and 6 columns
baseMean log2FoldChange lfcSE stat pvalue padj
<numeric> <numeric> <numeric> <numeric> <numeric> <numeric>
gene1 12.61910 -0.320613 0.728490 -0.440107 0.659860 0.992529
gene2 6.01265 -1.391883 0.926580 -1.502172 0.133053 0.972932
gene3 71.48063 -0.723134 0.319153 -2.265794 0.023464 0.805855
gene4 7.56355 -0.347392 0.715907 -0.485248 0.627501 0.992529
gene5 11.86764 0.275620 0.634890 0.434123 0.664199 0.992529
... ... ... ... ... ... ...
gene996 27.95109 -0.4044216 0.420641 -0.9614415 0.336330 0.992529
gene997 41.27173 -0.0273015 0.341122 -0.0800344 0.936210 0.996243
gene998 32.97030 0.1177037 0.357531 0.3292127 0.741995 0.992529
gene999 2.85416 -0.4862333 1.009720 -0.4815524 0.630124 0.992529
gene1000 3.24778 -0.2205322 1.003172 -0.2198348 0.826000 0.992529
> RESULT <- res[which(res$pvalue<0.05) & abs(res$log2FoldChange > 1.5),]
Warning message:
In which(res$pvalue < 0.05) & abs(res$log2FoldChange > 1.5) :
longer object length is not a multiple of shorter object length
> write.table(RESULT, file = "upregulated_test.txt")
> RESULT
log2 fold change (MLE): condition B vs A
Wald test p-value: condition B vs A
DataFrame with 30 rows and 6 columns
baseMean log2FoldChange lfcSE stat pvalue padj
<numeric> <numeric> <numeric> <numeric> <numeric> <numeric>
gene20 7.02072 1.88494 0.722398 2.60928 0.00907324 0.805855
gene45 4.76153 1.70103 0.877075 1.93943 0.05244909 0.851959
gene57 0.88648 3.34660 2.213826 1.51168 0.13061522 0.972932
gene74 1.27885 2.14873 1.632143 1.31651 0.18800332 0.992529
gene75 4.16297 1.79519 1.024415 1.75240 0.07970449 0.901962
... ... ... ... ... ... ...
gene868 1.34768 3.38540 1.50091 2.25557 0.0240976 0.805855
gene923 1.36849 1.68479 1.60413 1.05028 0.2935898 0.992529
gene932 6.63682 1.78308 0.86071 2.07164 0.0382986 0.851959
gene970 2.86467 1.79572 1.32983 1.35034 0.1769079 0.992529
gene977 4.02625 1.50587 1.01119 1.48920 0.1364334 0.974524
Iam also pasting the attaching the file which i created using above mentioned command.
Also needs to wrap
abs()
only around the LFC.