First, I’m not a statistician, so you may want to check with one. Also, I hope anyone reading this with more experience can correct or refine anything if/as necessary.

My main goal is to compare the DEGs identified in the different contrasts. For example, I want to know which genes are unique to each treatment and which genes are shared between treatments.

For this, you could start with set-operation analysis of per-contrast DEG lists: e.g., unions, intersections, set difference, mutual exclusivity (genes significant in only one contrast), complements, etc. However, be careful about interpreting results here. Across the planned contrasts, you’ll (almost certainly) have different effect sizes, variability, and thus statistical power. So, for example, if a DEG appears to be unique to one contrast, that doesn’t necessarily mean it is treatment-specific.

Could you explain what you meant by "it depends on the question you're trying to answer"? How would the answer differ if I am interested in DE genes for each contrast separately versus DE genes across all contrasts? How does that affect whether an additional multiple-testing correction across contrasts is needed?

When DESeq2 returns adjusted p-values, it’s going to be per contrast or coefficient (whatever you ran in results()).

For planned contrasts in one experiment, it’s fairly common to report the DESeq2 FDR-adjusted results per contrast without any additional across-contrast corrections. That would provide an answer to “how would the answer differ if I am interested in DE genes.” Those per-contrast DEG lists can then be compared with set operations.

What I meant by “it depends on the question” is that controlling FDR at 5% per contrast is not the same as controlling FDR at 5% across the whole set of contrasts.

If you wanted to treat all gene/contrast tests as one larger set of results, then a simple approach would be to pool the raw p-values from all planned comparisons and apply BH correction to the combined set. The same gene can appear more than once here; this is because testing the same gene in different contrasts represents different hypotheses. You would also need to be careful because DESeq2 does independent filtering separately for each results() call by default (i.e., it’s on by default but can be turned off), and padj can be set to NA for filtered genes.

Btw, I’m not recommending this. I’m just bringing it up because your initial post asked about whether an additional correction across comparisons would be needed. For your stated goal of comparing DEG lists, I’d probably start with the per-contrast DESeq2-adjusted results and then do the above set-operation analysis.

Actually, there are more formal hierarchical and stage-wise approaches for this kind of problem (e.g., pooling across contrasts and correcting). So, if you want or think it’s necessary to follow this route, I’d look into those. For example, see stageR.