I was trying to figure out how to perform Spearman rank test between expression value of a gene and a continuous outcome in DESeq2. One way I did was to extract the normalized count from counts(dds, normalized=TRUE), and then run cor.test in R. Alternatively, I put the continuous variable in the DESeq2 design formula and ran DESeq(dds). I could get the p and adjust p from the results.
It seems the latter approach could take advantage of the shrinkage estimator in DESeq2, and generate more significant genes. However, is there a way to get the Spearman correlation coefficient and p value from the DESeq2 output? Thanks!

We don't have a way to generate anything like a Spearman correlation in DESeq2.

You could alternatively use SAMseq, which has an implementation of Spearman designed specifically for RNA-seq count data.

See resp.type="Quantitative":

https://www.rdocumentation.org/packages/samr/versions/3.0/topics/SAMseq