degenerate-input-filtering

退化输入过滤指南

概述

退化输入是指不携带任何统计信息的数据点：恒定值特征、全为 NaN 的列、单序列比对、空文件及类似边缘情况。当这些数据进入统计检验或模型时，结果毫无意义——相关系数为 NaN、p 值为 1.0、得分为 0.0，或直接导致崩溃。本指南规定了强制做法：在任何分析之前检测并移除这类输入，并报告每项移除，以便下游使用者了解实际样本量。

关键概念

什么情况属于退化输入

如果数据点无法对正在计算的统计量做出贡献，就属于退化数据。其根本原因总是相同的：输入缺少该方法所需的变异性或完整性。

类型	示例	失败原因
恒定值特征	所有样本中表达量都相同的基因	方差 = 0；相关性、t 检验和倍数变化均无定义
全为 NaN 的特征	没有有效观测值的列	所有聚合结果均为 NaN
单序列比对	仅含一个序列的 BLAST 结果	得分 = 0.0；无法进行成对比较
空文件	0 字节的 FASTA 或 CSV 文件	解析器崩溃或返回空数据框
分组后仅有一个值	某个处理组中只有一个样本	组内方差无定义；组间比较毫无意义
零长度序列	空的 FASTA 条目	比对和 k-mer 工具会失败或产生无意义的结果
近乎恒定的特征	含一个离群值和 N-1 个相同值的基因	虽然方差不为零，但相关性会受到单个数据点的主导

为什么无提示的失败很危险

许多数值计算库不会因退化输入而报错，而是返回哨兵值：

• numpy.corrcoef 对恒定列返回 nan，且不发出警告。
• scipy.stats.spearmanr 在其中一个数组为常量时返回 (nan, nan)。
• scipy.stats.ttest_ind 在组的方差为零时返回 (nan, nan)。

这些 NaN 值会在管道中无声传播，污染聚合统计量、热图、火山图和基因排序列表。等研究人员发现问题时，可能已经无法确定是哪个上游步骤引入了 NaN。

报告义务

过滤却不报告几乎和不进行过滤一样有害。静默丢弃条目后，实际样本数或特征数会与用户预期不符，从而导致效能计算错误和汇总统计具有误导性。每个过滤步骤都必须打印：

1. 移除的条目数。
2. 移除原因。
3. 剩余的条目数。

决策框架

运行统计分析前，使用此决策树确定数据适用哪些过滤检查：

Is the input tabular (DataFrame)?
├── Yes
│   ├── Are there columns with zero unique values (all NaN)? → Remove, report count
│   ├── Are there columns with exactly one unique value (constant)? → Remove, report count
│   ├── Are there columns with fewer than N valid observations? → Remove, report count
│   └── Are there near-constant columns (1 outlier, rest identical)? → Flag or remove
└── No (file list, sequence set, etc.)
    ├── Are there empty files (0 bytes)? → Skip, report count
    ├── Are there single-entry collections (e.g., 1-sequence alignment)? → Skip, report count
    └── Are there zero-length entries within a file? → Filter entries, report count

场景	推荐检查	阈值
DE 分析前的基因表达矩阵	移除零方差基因，移除在少于 N 个样本中检测到的基因	N = max(3, 样本数的 10%)
两组特征之间的相关性分析	移除任一组中为常量的特征；要求至少有 >= 10 对共享的非 NaN 值	两个向量中的 nunique >= 2；共享观测值 >= 10
多序列比对评分	跳过序列数少于 2 的比对	序列数 >= 2
使用分组协变量的生存分析	移除事件数少于 2 的组	每组事件数 >= 2
对特征矩阵进行 PCA 或聚类	移除零方差和全为 NaN 的特征	方差 > 0，且至少有 1 个非 NaN 值
批次校正（例如 ComBat）	移除任何批次中缺失的特征	每个批次中均检测到该特征

最佳实践

1. 在进行任何统计计算之前过滤，而不是之后。 退化输入可能导致除以零、NaN 传播以及多重检验校正结果膨胀。测试完成后再过滤，损害已经造成。

2. 使用单个可复用的过滤函数。 集中处理逻辑可避免不同脚本之间出现不一致，并确保报告格式统一。请参阅下方工作流部分中的参考实现。

3. 在每个过滤步骤都打印计数。 即使计数为零，打印“移除了 0 项”也能确认检查已运行。这样日志便可自我说明并接受审计。

4. 设置符合领域特点的阈值，而不是依赖通用默认值。 基因表达分析可能要求至少在 10% 的样本中检测到，而缺失值较多的蛋白质组学数据集可能采用 5%。请记录阈值及其依据。

5. 谨慎处理近似常量特征。 某个特征只有一个值不同，严格来说方差不为零，但它与任何对象的相关性完全由这一个点决定。可考虑单独标记这些特征，而不是不加区分地将其纳入分析。

6. 切勿静默过滤。 删除数据却不报告，属于隐蔽的数据操纵。每次删除都必须记录原因和数量，确保有效样本量始终清晰可见。

7. 转换后重新检查。 对表达数据进行对数转换可能会因计数为零而产生负无穷值。合并数据集可能会引入新的 NaN 列。在任何转换或合并步骤后，都应再次进行退化检查。

常见陷阱

1. 对值恒定的列计算相关性，却没意识到结果是 NaN。 NaN 结果会悄无声息地传播到下游排序、热图或富集分析输入中，生成误导性图表。

   • 如何避免：计算相关性前，务必检查每列是否满足 nunique() >= 2。请使用下方的参考实现。

2. 对单序列比对进行评分，并将 0.0 分作为真实结果报告。 比对评分至少需要两个序列；只有一个序列的“比对”并非比对失败，而是根本没有进行比对。

   • 如何避免：评分前先统计序列数量。跳过并报告任何序列数少于 2 的比对文件。

3. 筛选数据却不报告数量，导致有效样本量令人困惑。 基于 15,000 个基因绘制的火山图与基于 18,500 个基因绘制的火山图会有很大差异，但如果没有日志消息，就没人知道实际使用的是哪种数据。

   • 如何避免：务必打印三项内容：移除数量、移除原因、剩余数量。每个笔记本和脚本都应包含这些信息。

4. 设置过低的最小样本数阈值，保留仅在 1-2 个样本中检测到的基因。 这类基因会导致方差估计不稳定、p 值不可靠，并抬高假发现率。

   • 如何避免：将 max(3, int(0.1 * n_samples)) 设为下限。对于样本量较小的研究，如果 10% 仍少于 3，则应调高阈值。

5. 合并或转换数据后忘记再次检查。 合并两个数据集可能会因特征不重叠而产生全为 NaN 的列。对零值进行对数转换会产生负无穷值。这两种情况都属于退化数据。

   • 如何避免：每次合并、连接或数学转换后，都要再次运行退化数据筛选。

6. 使用 df.var() > 0，却没有正确处理 NaN。 如果某列全为 NaN，var() 会返回 NaN，而 NaN 不满足 > 0，因此该列会被删除——但日志记录的原因是“零方差”而非“全为 NaN”，这会造成误导。

   • 如何避免：检查零方差之前，先单独检查全为 NaN 的列。分别使用各自的标签报告每种类别。

7. 因为数据来自经过整理的数据库，就认为输入是干净的。 公共数据库中包含占位值、缺失条目和已撤回记录。即使来源可信，也务必进行验证。

   • 如何避免：将所有输入都视为不可信数据。无论数据来源如何，都要执行完整的退化数据检查。

工作流

在进行任何统计分析之前，都应按顺序执行以下流程。

1. 步骤 1：加载并检查

   • 将数据加载到 DataFrame 或等效结构中。
   • 打印数据形状：行数、列数和 dtype 概览。

2. 步骤 2：移除全为 NaN 的特征

   • 找出所有值均为 NaN 的列（或行，具体取决于数据方向）。
   • 移除这些列或行，并打印数量。

3. 步骤 3：移除值恒定的特征

   • 找出 nunique() <= 1 的列（不计 NaN）。
   • 移除这些列，并打印数量。

4. 步骤 4：移除有效观测值过少的特征

   • 设置阈值（例如，至少占样本数的 10% 或至少为 3）。
   • 移除低于阈值的特征，并打印数量。

5. 步骤 5：应用领域特定检查

   • 对序列数据：跳过单序列比对、空文件和零长度条目。
   • 对表达数据：筛选低检出率基因。
   • 对相关性输入：验证两个向量都不是常量，并且具有足够数量的共同观测值。

6. 步骤 6：报告摘要

   • 打印最终摘要：输入总数、按类别统计的移除总数、剩余总数。

参考实现

import pandas as pd
import numpy as np

def filter_degenerate(data, context="features"):
    """Filter degenerate data points and report what was removed.

    Always call this BEFORE statistical analysis.
    """
    n_before = len(data)

    # Filter based on data type
    if isinstance(data, pd.DataFrame):
        # Remove constant columns (zero variance)
        non_constant = data.loc[:, data.nunique() > 1]
        # Remove all-NaN columns
        non_nan = non_constant.dropna(axis=1, how='all')
        filtered = non_nan
    elif isinstance(data, list):
        # Remove None, empty, and zero-length items
        filtered = [x for x in data if x is not None and len(x) > 0]
    else:
        filtered = data

    n_after = len(filtered) if not isinstance(filtered, pd.DataFrame) else filtered.shape[1]
    n_removed = n_before if not isinstance(data, pd.DataFrame) else data.shape[1]
    n_removed = n_removed - n_after

    # MANDATORY: Print the count
    print(f"Degenerate {context} filtered: {n_removed} / {n_before} removed")
    print(f"Remaining {context}: {n_after}")

    return filtered

序列比对筛选

# Filter single-sequence alignments (score=0.0 is meaningless)
alignment_scores = {}
for aln_file in alignment_files:
    n_seqs = count_sequences(aln_file)
    if n_seqs < 2:
        print(f"Skipping {aln_file}: only {n_seqs} sequence(s)")
        continue
    score = compute_alignment_score(aln_file)
    alignment_scores[aln_file] = score

print(f"Filtered {len(alignment_files) - len(alignment_scores)} single-sequence alignments")

基因表达筛选

# Filter genes with zero variance (constant expression)
gene_vars = expression_df.var(axis=1)
zero_var = (gene_vars == 0).sum()
print(f"Genes with zero variance: {zero_var}")
expression_filtered = expression_df.loc[gene_vars > 0]

# Filter genes detected in too few samples
min_samples = max(3, int(0.1 * expression_df.shape[1]))  # At least 10% of samples
detected = (expression_df > 0).sum(axis=1)
low_detection = (detected < min_samples).sum()
print(f"Genes detected in < {min_samples} samples: {low_detection}")
expression_filtered = expression_filtered.loc[detected >= min_samples]

相关性预筛选

from scipy.stats import spearmanr

# Filter features before computing correlations
valid_correlations = {}
skipped_constant = 0
skipped_few_values = 0

for feature in features:
    x = data_x[feature].dropna()
    y = data_y[feature].dropna()

    # Check for constant values
    if x.nunique() < 2 or y.nunique() < 2:
        skipped_constant += 1
        continue

    # Check for sufficient data points
    common = x.index.intersection(y.index)
    if len(common) < 10:
        skipped_few_values += 1
        continue

    rho, pval = spearmanr(x[common], y[common])
    valid_correlations[feature] = rho

print(f"Skipped (constant values): {skipped_constant}")
print(f"Skipped (< 10 valid pairs): {skipped_few_values}")
print(f"Valid correlations computed: {len(valid_correlations)}")

预期输出格式

格式良好的筛选报告如下所示：

Input: 18500 genes
Filtered: 230 genes with zero variance
Filtered: 45 genes with < 3 valid samples
Remaining: 18225 genes for analysis

延伸阅读

• scipy.stats 文档 —— 统计函数及其处理边缘输入的行为参考，包括常量数组和 NaN 值。
• pandas DataFrame.dropna —— 关于通过 axis 和阈值控制移除 NaN 的官方文档。
• numpy NaN 处理 —— NumPy 文档，介绍算术和聚合函数中 NaN 的传播语义。
• Lun, A. 等（2016）《使用 Bioconductor 对单细胞 RNA-seq 数据进行低层级分析的逐步工作流》 —— 展示基因层面的质量控制过滤，包括在归一化和差异表达分析前移除低丰度和零方差基因。
• Law, C.W. 等（2018）《使用 limma、Glimma 和 edgeR 进行 RNA-seq 分析轻而易举》 —— 介绍类似 filterByExpr 的基因过滤方法，用于在差异表达分析前移除低表达和缺乏信息的基因。

相关技能

• nan-safe-correlation —— 处理缺失值时计算相关性的技巧；可补充本指南中对相关性输入进行预过滤的步骤。
• statistical-analysis —— 关于选择和应用统计检验的更广泛指导；假设已按照本指南移除退化输入。