一、公式声明
需要声明，这一实现使用的标准组的Ridit方差为贝塞尔校正版本，而卡方统计量的公式采用如下形式：
——式子1
其中：
注意以下公式默认了有序评分对应的隐连续得分是服从均匀分布的，因此才将，然后把式子1化简成了以下形式：
——式子2
但这一假设通常不成立，因此几乎总是需要有结校正。有结校正本质上就是将式子2转换回了式子1。因此我们总是采用式子1计算卡方统计量。
二、R代码实现
ridit.test = function( data, std_group_ratio = 3, conf = 0.95 ) { # 使用双侧检验，比较不同组请看图 if (is.null(nrow(data))) { stop("检验数据中至少包含两个对象！") } sample_size <- sum(data) group_sample_sum <- rowSums(data) std_group_idx <- -1 df <- nrow(data) - 1 group_names <- if (!is.null(rownames(data))) { rownames(data) } else { paste("G", 1:nrow(data)) } # 1. 确定标准组 sort_sample_sum <- sort(group_sample_sum, decreasing = TRUE) if (sort_sample_sum[1] >= 3 * sort_sample_sum[2]) { std_group_idx <- which.max(group_sample_sum) group_sample_sum <- group_sample_sum[-std_group_idx] std_group <- data[std_group_idx, ] control_group <- data[-std_group_idx, ] std_group_name <- group_names[std_group_idx] control_names <- group_names[-std_group_idx] } else { std_group <- colSums(data) control_group <- data std_group_name <- "merged" control_names <- group_names } # 2. 计算每个评分映射到的Ridit值 N_std <- sum(std_group) cum_freq <- c(0, cumsum(std_group)[-length(std_group)]) std_ridit <- (cum_freq + std_group / 2) / N_std # 3. 计算对照组的平均Ridit值 control_avg_ridit <- (control_group / group_sample_sum) %*% std_ridit # 4. 计算标准组的方差 std_variance <- ((1 / N_std) * sum(std_group * std_ridit^2) - 1/4) * (N_std / (N_std - 1)) control_SE <- sqrt(std_variance / group_sample_sum) # 5. 根据对照组数量的不同安排不同检验 if (is.null(nrow(control_group)) || nrow(control_group) == 2) { if (is.null(nrow(control_group))) { Z <- (control_avg_ridit - 0.5) / control_SE } else { avg_R_diff <- abs(control_avg_ridit[1] - control_avg_ridit[2]) combine_SE <- sqrt(sum(1 / group_sample_sum) * std_variance) Z <- avg_R_diff / combine_SE } statistic <- c(Z = Z) p.value <- pnorm(Z, lower.tail = FALSE) * 2 } else { Z_list <- (control_avg_ridit - 0.5) / control_SE W <- sum(Z_list^2) statistic <- c(W = W) p.value <- pchisq(W, df, lower.tail = FALSE) } # 6. 画图 Z_line <- qnorm(1 - (1 - conf) / 2) control_ci_lower <- control_avg_ridit - Z_line * control_SE control_ci_upper <- control_avg_ridit + Z_line * control_SE control_groups_num <- max(nrow(control_group), 1) margin <- 0.2 if (control_groups_num == 1) { x_positions <- 1.5 xlim_left <- 0.5 xlim_right <- 2.5 } else { # 计算位置，让线条均匀分布但不在最边缘 total_width <- control_groups_num - 1 start <- 1 + margin end <- control_groups_num - margin x_positions <- seq(start, end, length.out = control_groups_num) xlim_left <- 0.5 xlim_right <- control_groups_num + 0.5 } plot( 0, 0, ylim = c(0, 1), xlim = c(xlim_left, xlim_right), xlab = "group", ylab = "R score", main = "Ridit value confidence interval", col = "gray7", xaxt = "n", yaxt = "n" ) abline(h = 0.5, col = "black", lty = 2, lwd = 1.5) for (i in 1:control_groups_num) { lines( c(i, i), c(control_ci_lower[i], control_ci_upper[i]), col = "black", lwd = 2 ) } # 添加X轴标签 axis(1, at = 1:control_groups_num, labels = control_names) axis(2, at = 0.5, labels = std_group_name) result <- list( statistic = statistic, df = df, p.value = p.value, parameter = c(df = df), data.name = deparse(substitute(data)), std_variance = std_variance, std_redit = std_ridit, control_avg_ridit = control_avg_ridit, method = "ridit Test", control_group = control_group, std_group = std_group ) class(result) <- "htest" return(result) }三、测试结果
3.1单样本测试
test1 = matrix( c( 99, 193, 65, 28, 12, 14, 45, 17, 20, 7 ), nrow = 2, byrow = TRUE) result <- ridit.test(test1) result输出结果：
ridit Test data: test1 Z = 4.3568, df = 1, p-value = 1.319e-05详细信息：
> result$std_variance [1] 0.0722512 > result$std_redit [1] 0.1246851 0.4924433 0.8173804 0.9345088 0.9848866 > result$control_avg_ridit [,1] [1,] 0.6153921因此拒绝原假设，认为G1和G2存在显著性差异。且根据图示结果，G1的得分要大于G2，注意解释时要看累积概率即Ridit值的计算方向。
3.2 双样本测试
test2 = matrix( c( 6, 18, 19, 27, 25, 15, 31, 31, 32, 19 ), nrow = 2, byrow = TRUE) result <- ridit.test(test2) result输出结果：
ridit Test data: test2 Z = 2.488, df = 1, p-value = 0.01285因为组1和组2样本差距不大，因此将合并组的结果作为标准组，即纵轴的merged；检验结果表明在组1和组2在的情况下存在显著性差异；从图中可以看出两个组的置信区间几乎就要有重叠部分。
3.3 多样本测试
test3 = matrix( c( 15, 45, 153, 231, 56, 39, 89, 198, 126, 48, 89, 177, 134, 88, 12, 13, 23, 86, 257, 121 ), nrow = 4, byrow = TRUE) result <- ridit.test(test3) result输出结果：
ridit Test data: test3 W = 441.51, df = 3, p-value < 2.2e-16检验结果表明4个组存在非常显著的差异。从图中可以看出，合并组作为了标准组，四个组中，组1和组4的Ridit得分要显著高于组2和组3；
)
