1. Support
2. Help - Theory & Case Study
3. How Can We Find Important Genes That Do Not Meet DEG Criteria in RNA-Seq?

Is Current RNA-Seq Analysis Too Focused on Downstream Changes?

In RNA-Seq analysis, differential expression analysis using DESeq2 or edgeR is widely used, often followed by enrichment analyses such as GO analysis and pathway analysis. These approaches are highly useful for identifying genes that show large expression changes between conditions and for examining biological functions shared by those genes. However, there is growing concern that looking only at genes with large expression changes may make it difficult to identify genes that lie closer to the upstream causes of those changes.

A 2026 review published in Briefings in Bioinformatics, Multimodal bioinformatic analyses of genome-scale expression beyond gene-centric differential expression , summarizes the development of approaches that go beyond conventional gene-centric differential expression analysis by incorporating gene co-expression networks, gene regulatory networks, multi-omics integration, machine learning, and other methods to understand gene expression in a broader regulatory context. The review describes how analysis is expanding from single-gene differential expression toward systems-level approaches that combine network information with multiple omics modalities.

This may be a somewhat rough way to put it, but one could say that the field is increasingly moving toward looking not only at DEGs, but also at genes that are not selected as DEGs (non-DEGs), as well as genes whose expression changes in a direction opposite to what might have been expected.

Finding Genes Closer to the Cause Is Not So Simple

Even when we try to identify genes that are closer to the cause of a biological phenomenon of interest, this is far from straightforward in practice.

In a previous case study, we used ChIP-Atlas to search for candidate upstream regulators by identifying transcription factors with binding peaks near genes that changed in RNA-Seq. However, the candidate transcription factors rarely showed the expression changes we would have expected, and the candidates as a group did not even consistently shift toward either increased or decreased expression. This illustrates that, although methods for inferring genes closer to the causal mechanism are needed, doing so remains a very difficult challenge in practice.

A more realistic approach may be to combine RNA-Seq with other types of experimental data through multi-omics integration. In the previous two case studies, we examined studies that combined RNA-Seq and RIP-seq. RIP-seq can be used to identify RNAs associated with a specific RNA-binding protein (RBP) by immunoprecipitating protein-bound RNAs and sequencing them.

What was particularly interesting was that, for RNA groups defined by RIP-seq, even when the expression changes of individual genes were too small to exceed conventional DEG thresholds, the distribution of the gene group could still show a statistically significant shift relative to a control group in a KS test.

For details, see “Some Genes That Were Not Significant in DEG Analysis Have Been Functionally Validated Experimentally” and “Have RNA-Seq Analysis Methods Been Too Biased? Exploring Possibilities Beyond DEG Analysis”.

This Time, We Combine Multi-omics with Network Analysis

This time, we would like to try a somewhat different approach.

As in GSE124280, which we examined in the previous two case studies, we will use another dataset in which both RNA-Seq and RIP-seq were obtained: GSE86590.

This time, however, we will try an approach that

combines experimental multi-omics data with algorithm-based network analysis.

The network analysis method used here is DENetwork . DENetwork was published in Nucleic Acids Research in 2025. It is designed to consider not only DE genes but also non-DE genes that may be overlooked by conventional differential expression analysis, and to evaluate candidate genes based on their contribution to information flow within a gene network. One of the motivations described in the paper is that DE methods such as DESeq2 focus primarily on changes in expression level and may therefore miss genes that show only small expression changes but nevertheless play important roles in a network.

However, we do not regard DENetwork as a method for “identifying the causal gene.” What we want to examine here is what happens when genes selected by DENetwork from a PPI network as candidates bridging the experimental perturbation of YBX2 knockout and the resulting differential expression changes are compared with RIP-seq data obtained through a different experimental principle.

Dataset: YBX2 Knockout Brown Adipose Tissue

GSE86590 investigates the role of the RNA-binding protein YBX2 in brown adipose tissue (BAT). BAT RNA-Seq data are available from wild-type and Ybx2 knockout mice at room temperature and after cold exposure. YBX2 RIP-seq and an IgG control were also obtained to examine RNAs associated with YBX2.

In the original study, RNA-Seq and RNA immunoprecipitation were integrated to identify YBX2 mRNA targets, including Pgc1α, and YBX2 was shown to be involved in the regulation of RNA stability during cold-induced BAT activation.

In other words, this dataset contains two types of data that are well suited to our purpose:

• RNA-Seq after loss of YBX2
• RIP-seq measuring RNAs associated with YBX2

Connecting YBX2 and Down-DEGs with DENetwork

We used YBX2 alone as the source and genes downregulated in Ybx2 knockout mice after cold exposure as the target Down-DEGs. In other words, we searched for paths connecting YBX2 to Down-DEGs in a PPI network and evaluated the internal nodes located between them.

With K=3, which allows a maximum of three nodes including the source and target, or at most one internal node, no qualifying path was obtained.

With K=4, which allows up to two internal nodes, two paths were obtained:

• YBX2 → PCBP1 → KCTD13 → FASN
• YBX2 → PCBP1 → TMEM173 → FASN

With K=6, which allows up to four internal nodes, the following five paths remained. Interestingly, all of the final paths consisted of four nodes, with two internal nodes between the source and target.

1. YBX2 → PCBP1 → TMEM173 → FASN
2. YBX2 → PCBP1 → KCTD13 → FASN
3. YBX2 → PCBP1 → PPP1R9B → FASN
4. YBX2 → PCBP1 → GRN → FASN
5. YBX2 → PCBP1 → STAU1 → FASN

The six internal nodes were PCBP1, TMEM173, KCTD13, PPP1R9B, GRN, and STAU1. PCBP1 was shared by all of the paths. If we were to investigate this result further, one of the first things to check would be whether PCBP1 is a hub in the PPI network used here. Genes with many connections in a PPI network may appear in many paths simply because of their high connectivity, potentially introducing bias into network-based analysis.

We next examined the expression changes of these genes. In the YBX2 knockout samples, Ybx2 itself and Fasn showed large decreases in expression, to approximately 0.13-fold and 0.04-fold, respectively. In contrast, the internal nodes showed the following changes:

PCBP1 decreased slightly, whereas KCTD13, PPP1R9B, GRN, and STAU1 increased slightly. TMEM173 increased by approximately 2.8-fold. In other words, if we had followed only the genes that decreased after YBX2 knockout, none of these six genes would have been selected as candidates. DENetwork instead identified genes that would not have been selected based solely on the direction or magnitude of their expression changes, placing them as internal nodes between YBX2 and the Down-DEGs.

However, this should not be interpreted as proof of a biological pathway such as YBX2 → PCBP1 → TMEM173 → FASN. The PPI network used by DENetwork is based on previously reported protein-protein interactions, and the presence of such a path does not demonstrate that the same sequence of interactions actually functions as a biological pathway in BAT under the conditions examined here.

Evaluating the Six DENetwork Candidates with RIP-seq

We compared YBX2 RIP-seq with IgG RIP-seq as a control. All six DENetwork internal nodes satisfied the criterion of YBX2 RIP FPKM >1 and a YBX2/IgG ratio of at least 2-fold.

Under more stringent ratio thresholds:

• YBX2/IgG ratio of at least 2-fold: 6/6
• YBX2/IgG ratio of at least 4-fold: 5/6
• YBX2/IgG ratio of at least 8-fold: 3/6

For TMEM173, the IgG RIP FPKM was 0, so the YBX2/IgG ratio itself could not be calculated. However, the YBX2 RIP FPKM was 1.53, so for the purpose of this threshold-based classification, we included TMEM173 in the “8-fold or higher” group.

At first glance, this may appear to show substantial overlap between the genes selected by DENetwork and the YBX2 RIP-seq results. However, the background frequency also needs to be considered. We therefore used the 8,576 genes shared between the PPI network and the RNA-Seq dataset that were eligible for the DENetwork analysis as the overall analysis set, and compared the six internal nodes with the remaining genes.

When RIP FPKM >1 was used as the signal threshold:

A relatively high proportion of the internal nodes showed RIP enrichment. However, YBX2 RIP enrichment was also common among the background genes. Therefore, the enrichment among these six internal nodes was not statistically significant.

When a stronger RIP signal threshold of YBX2 RIP FPKM >5 was used:

With RIP FPKM >5 and a YBX2/IgG ratio of at least 2-fold, the p value was approximately 0.049. However, there were only six internal nodes. In addition, we compared multiple cutoff conditions, and this p value was not adjusted for multiple comparisons. This result should therefore not be treated as strong statistical evidence. Furthermore, there is no universal RIP-seq cutoff defining what FPKM value should be regarded as a sufficiently strong RIP signal. For this reason, we show results using both RIP FPKM >1 and RIP FPKM >5 to examine how much the results depend on the cutoff.

Intuitively, one might expect that the more stringent the YBX2/IgG ratio threshold becomes—from 2-fold to 4-fold to 8-fold—the more strongly the DENetwork internal nodes would be enriched. However, the opposite pattern was observed. For example, with RIP FPKM >1, increasing the threshold from 2-fold to 8-fold reduced the proportion of background genes satisfying the criterion from 65.3% to 35.1%. At the same time, the number of internal nodes satisfying the criterion also decreased from 6/6 to 3/6. In other words, the stricter cutoff removed internal-node genes as well as background genes, so the statistical difference became weaker rather than stronger.

These results do not indicate that DENetwork preferentially selects RNAs showing the strongest RIP enrichment. However, this is not necessarily surprising. RIP-seq examines associations between YBX2 protein and RNA, whereas the PPI network used by DENetwork represents interactions between proteins. These are fundamentally different biological relationships, so the strength of RIP enrichment would not necessarily be expected to correspond directly to the DENetwork results.

What may be more important is whether two different types of evidence converge on the same candidate genes, and whether that overlap can help narrow down genes that are functionally important. In this article, we cannot perform biological functional experiments on these candidate genes, so we cannot determine whether they are truly functionally important. Nevertheless, combining these approaches may be worth considering as one way to search for candidates that are difficult to identify through DEG analysis alone.

It Is Not a Choice Between Algorithms and Multi-omics

When trying to identify genes that are closer to the causes of biological changes in RNA-Seq analysis, the discussion can easily become focused on the question, “Which algorithm should we use?” However, we think it may be more realistic to explore candidates by combining multi-omics data with analytical methods than to try to narrow down causal candidates using RNA-Seq analysis alone.

Of course, there are limitations. In this case study, for example, questions remain about the appropriate cutoff for RIP-seq, interactions in the PPI network that may not match the experimental conditions, bias toward hub genes, and dependence on the models and parameters used by the analysis algorithm. In other words, neither a single dataset nor a single algorithm is definitive.

That is why there may be value in asking whether different types of evidence point to the same candidates. How those candidates should be interpreted and what should be tested next may also benefit from open discussion among researchers with different areas of expertise.

The direction described in the 2026 review as moving “beyond gene-centric differential expression” does not necessarily mean abandoning DEG analysis. Rather, it may be more appropriate to view DEG analysis as a starting point and then integrate additional information from networks, multi-omics, machine learning, and other approaches as RNA-Seq analysis continues to expand.