基因芯片数据分析-2: Quantile normalizati

分位数归一化就是将多种分布映射到同一种分布。可以指定一个参考分布，将其映射到其它的统计分布；或者从需要处理的这些不同数据中计算出一个参考分布。

原理实现：

1. 存在一个原始数组如下，每一列代表一个分布，最终我们将它们转换为同一个分布

A    5    4    3
B    2    1    4
C    3    4    6
D    4    2    8

2. 原始数组每列中按从小到大排序

A 2 1 3
B 3 2 4
C 4 4 6
D 5 4 8

3. 获得原始数组中每个元素在所属列中按从小到大的排序之后的序号

A    iv    iii   i
B    i     i     ii
C    ii    iii   iii
D    iii   ii    iv

4. 获取原始数组每列中按从小到大排序之后每行的算术平均数

A (2 1 3)/3 = 2.00 = rank i
B (3 2 4)/3 = 3.00 = rank ii
C (4 4 6)/3 = 4.67 = rank iii
D (5 4 8)/3 = 5.67 = rank iv

5. 利用第4步的对应关系，对第三部的数据进行替换

A    5.67    5.17    2.00
B    2.00    2.00    3.00
C    3.00    5.17    4.67
D    4.67    3.00    5.67

6. 查看每列数据的统计情况

Min.   :2.000   Min.   :2.000   Min.   :2.000  
1st Qu.:2.750   1st Qu.:2.750   1st Qu.:2.750  
Median :3.833   Median :4.083   Median :3.833  
Mean   :3.833   Mean   :3.833   Mean   :3.833  
3rd Qu.:4.917   3rd Qu.:5.167   3rd Qu.:4.917  
Max.   :5.667   Max.   :5.167   Max.   :5.667  

分位数归一化这种方法常用于芯片测序数据，但这种方法也过于粗暴，结果往往不是很好；在芯片数据预处理时，可以使用R中preprocessCore包中打包好的分位数归一化函数。

# 安装

if (!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("preprocessCore")

# 使用

library("preprocessCore")
data <- data.frame(a = sample(c(1:20),5),
                b =  sample(c(1:20),5),
                c = sample(c(1:20),5)
)
data
   a  b  c
1  7 10 11
2 20 11 16
3  6  8  7
4  8 13  9
5 17  3 14
normalize.quantiles(x=as.matrix(df))
          [,1]      [,2]      [,3]
[1,]  8.000000  9.666667  9.666667
[2,] 16.333333 14.000000 16.333333
[3,]  5.333333  8.000000  5.333333
[4,]  9.666667 16.333333  8.000000
[5,] 14.000000  5.333333 14.000000

normalize.quantiles()使用的参考分布是数据的算术平均数；需要指定参考分布，可以参考DAVE TANG'S BLOG: Quantile normalisation in R从头写的代码

下面是我稍微改动了一下的代码，添加指定指定参考分布的参数。

df <- data.frame(one=c(5,2,3,4),
                 two=c(4,1,4,4),
                 three=c(3,4,6,8)
)
quantile_normalisation <- function(df, ref_distribute){
  df_rank <- apply(df,2,rank,ties.method="min")
  df_sorted <- data.frame(apply(df, 2, sort))
  
  if(missing(ref_distribute)){
    df_mean <- apply(df_sorted, 1, mean)
  }else{
    df_mean <- ref_distribute
  }
  
  index_to_mean <- function(my_index, my_mean){
    return(my_mean[my_index])
  }
  
  df_final <- apply(df_rank, 2, index_to_mean, my_mean=df_mean)
  rownames(df_final) <- rownames(df)
  return(df_final)
}
quantile_normalisation(df)
       one      two    three
1 5.666667 3.666667 2.000000
2 2.000000 2.000000 3.666667
3 3.666667 3.666667 4.666667
4 4.666667 3.666667 5.666667

test = c(1,2,3,4)
quantile_normalisation(df, ref_distribute=test)
  one two three
1   4   2     1
2   1   1     2
3   2   2     3
4   3   2     4

参考:
Quantile normalization
DAVE TANG'S BLOG: Quantile normalisation in R