代谢组学分析中PQN归一化代码运行全为NA的解决咨询

我的数据前6行如下：

我希望使用PQN方法进行归一化，以下是使用的R代码：

library(dplyr)

x <- read.csv("mydata.csv", header=TRUE)

normalize_data_pqn <-
  function(x,
           pqn_reference = "median",
           pgn_reference_sample = NULL) {
    
    if (!is.null(pgn_reference_sample)) {
      if (pqn_reference == "mean") {
        mx <- as.numeric(apply(x[, pgn_reference_sample], 1 , mean))
      }
      if (pqn_reference == "median") {
        mx <- as.numeric(apply(x[, pgn_reference_sample], 1 , median))
      }
    } else {
      if (pqn_reference == "mean") {
        mx <- as.numeric(apply(x, 1 , mean))
      }
      if (pqn_reference == "median") {
        mx <- as.numeric(apply(x, 1 , median))
      }
    }
    
    # do the actual normalization
    new_x <-
      purrr::map(
        x,
        .f = function(value) {
          as.numeric(value / median(as.numeric(value / mx)))
        }
      ) %>%
      dplyr::bind_cols()
    
    rownames(new_x) <- rownames(x)
    return(new_x)
  }

normalize_data_pqn(x, pqn_reference = "median", pgn_reference_sample = NULL)

运行上述代码后，所有结果值变为NA，且出现警告信息：

问题原因及解决办法 #

1. 数据格式与非数值列干扰 #

通常代谢组数据的第一列是特征名称（字符型），但原代码直接读取整个数据集传入函数，非数值列会导致计算时生成NA。

• 解决：读取时指定行名，或提前筛选数值列：

  # 方法1：读取时将第一列设为行名
  x <- read.csv("mydata.csv", header=TRUE, row.names=1)
  # 方法2：筛选所有数值型列
  x_num <- dplyr::select_if(x, is.numeric)
  

2. 0值导致的计算错误 #

数据中存在0值时，value/mx会生成Inf或NaN，中位数计算失败最终返回NA。

• 解决：给所有数值加极小值避免0，同时过滤掉参考值mx为0的特征：

  # 给数值加极小值
  x_num <- x_num %>% mutate_all(~ . + 1e-10)
  

3. 依赖包缺失 #

原函数使用purrr包但未加载，会导致map函数报错。

• 解决：开头加载purrr：

  library(dplyr)
  library(purrr)
  

4. 函数逻辑补全 #

原函数未处理na.rm，且未过滤无效特征，补充后可避免NA扩散：

修正后的完整代码 #

library(dplyr)
library(purrr)

# 读取并预处理数据
x <- read.csv("mydata.csv", header=TRUE, row.names=1)
x_num <- select_if(x, is.numeric) %>% mutate_all(~ . + 1e-10)

normalize_data_pqn <-
  function(x,
           pqn_reference = "median",
           pgn_reference_sample = NULL) {
    
    # 计算特征参考值（每行的中位数/均值）
    if (!is.null(pgn_reference_sample)) {
      if (pqn_reference == "mean") {
        mx <- as.numeric(apply(x[, pgn_reference_sample], 1 , mean))
      } else if (pqn_reference == "median") {
        mx <- as.numeric(apply(x[, pgn_reference_sample], 1 , median))
      }
    } else {
      if (pqn_reference == "mean") {
        mx <- as.numeric(apply(x, 1 , mean))
      } else if (pqn_reference == "median") {
        mx <- as.numeric(apply(x, 1 , median))
      }
    }
    
    # 过滤参考值为0的无效特征
    valid_idx <- mx != 0
    x_valid <- x[valid_idx, ]
    mx_valid <- mx[valid_idx]
    
    # 执行PQN归一化，忽略NA值
    new_x <-
      purrr::map(
        x_valid,
        .f = function(value) {
          ratio <- value / mx_valid
          scale_factor <- median(ratio, na.rm = TRUE)
          as.numeric(value / scale_factor)
        }
      ) %>%
      dplyr::bind_cols()
    
    # 补回无效特征的NA行（如果存在）
    if (sum(!valid_idx) > 0) {
      na_rows <- matrix(NA, nrow = sum(!valid_idx), ncol = ncol(x)) %>% 
        as.data.frame() %>% 
        setNames(colnames(x))
      new_x <- rbind(na_rows, new_x)
    }
    
    rownames(new_x) <- rownames(x)
    return(new_x)
  }

# 运行归一化
normalized_data <- normalize_data_pqn(x_num, pqn_reference = "median", pgn_reference_sample = NULL)

内容的提问来源于stack exchange，提问作者Fang