3.4.10. 中位数问题

中位数(median),又称中值,是一个重要的统计指标。它是将一组数据按顺序排列后居于中间位置的数,可将数据划分为个数相等的两部分,其中一个部分中所有数都大于等于中位数,另一个部分中所有数都小于等于中位数。对于有限的数集,如果其中有奇数个数,那么把所有数排序后正中间的那一个数的值就是中位数;如果有偶数个数,排序后最中间两个数的平均值是中位数。

示例

数列 [7,1,3,20,5] 的中位数是5。数列中比5小的数和比5大的数个数相等,都是2个。如果把这5个数排序后得到 [1,3,5,7,20],可以发现5位于正中间的位置。

数列 [1,3,10,2] 的中位数是2和3的平均数2.5。数列中比2.5小的数和比2.5大的数都是2个。如果把数列排序为 [1,2,3,10],就可以直观地看出2和3是最中间的两个数,所以中位数取它们的平均数。

简而言之,一个数列的中位数就是将其排序后,位于最中间位置的数值。对于这个数列而言,其中比中位数小的数和比中位数大的数应当是数量相等的。所以说中位数可以将数列中的数据划分为个数相等的两个部分,一个部分中所有数都大于等于中位数,另一个部分中所有数都小于等于中位数。

关于中位数,有一点概念必须澄清:它是一个单独的数,而不是数列中的一个元素!

如果给出一组数值,求它的中位数是很简单的,只需要对数组进行排序,然后根据元素个数的奇偶性来进行计算即可。

C++提供了一个非常方便实用的快速排序函数 std::sort(),可以用来对数组中指定的一段元素进行排序。调用 sort() 函数,只需要提供数组中要排序的那一段元素的头尾指针即可。要使用这个函数,需引入 algorithm 库并使用命名空间 std,下面是用法示例:

#include <algorithm>
using namespace std;

int a[100];
for (int i = 0; i < 100; ++i) scanf("%d", &a[i]); // 读入100个元素

sort(a, a + 10);      // 对 a[0] 到 a[9] 之间的10个元素进行排序,注意含头不含尾的范围表示原则
sort(a + 20, a + 40); // 对 a[20] 到 a[39] 之间的20个元素进行排序
sort(a, a + 100);     // 对 a[0] 到 a[99],即整个数组进行排序

所以,假如有一个题目,要读入 \(n\le1000\) 个整数并求其中位数,那么程序将非常简单:

#include <cstdio>
#include <algorithm>

using namespace std;

int main()
{
        int n, a[1001];
        scanf("%d", &n);
        for (int i = 0; i < n; ++i) scanf("%d", &a[i]);  // 读数
        sort(a, a + n);  // 排序
        // 根据 n 的奇偶性来计算中位数并输出,保留2位小数
        printf("%.2lf\n", n % 2 ? (double) a[n/2] : (a[n/2] + a[n/2-1]) / 2.0);

}

这个程序看懂了吗?看懂了我们就继续讨论一个更加复杂的问题。

现在有两个有序数组 \(A\)\(B\),各有 \(m\)\(n\) 个元素,两个数组都已经完成升序排序,即元素从小到大排列。现在要求快速求出这两个数组中所有元素的中位数,时间复杂度不高于 \(O\bigl(\log(m+n)\bigr)\)。其中 \(0\le m,n \le 10^6\),保证两个数组不会同时为空。

3.4.10.1. 归并法求两个有序数组的中位数

如果两个数组都是乱序的,为了求出它们所有元素的中位数,我们可以把两个数组的元素都合到一个大数组里,然后对大数组进行排序来求出中位数。但排序算法最快的时间复杂度为 \(O\bigl((n+m)\log(n+m)\bigr)\),远大于问题所要求的上限 \(O\bigl(\log(n+m)\bigr)\)。这个方法的主要问题出在没有利用两个原始数组都已经有序这个条件。

那么即然两个数组都已经排好序了,我们可以用一种叫做归并的操作把两个数组中的元素保持有序地合到另一个大数组里去,时间复杂度为线性的 \(O(n+m)\)。归并操作非常重要,在大量的算法问题中会用到,基本排序算法中最重要的归并排序也是基于归并来完成的,一定要熟练掌握。

归并的思路来源于整队。比如有两列队伍,都是按照身高从矮到高排列好的,现在要把两列队伍合并为一列,并且保持身高从矮到高排列。这个整队过程非常简单,从两列队伍各自排在队头的人开始,每次挑出两个队头中比较矮的那个人出列,排进新的队伍的尾部,如果队头两个人身高一样就任选其中一人即可。不断循环这一过程,直到某个队伍的人全部挑完,把另一个队伍中剩余的人按原顺序排到新队伍的后面就完成了。下面举个简单的例子来说明这个过程:

初始状态:

队伍一: [1, 3]

队伍二: [2, 6, 7]

新队伍: []

第1步:队伍一的队头1较小,挑出来放到新队伍尾部

队伍一: 1, [3]

队伍二: [2, 6, 7]

新队伍: [1]

第2步:队伍二的队头2较小,挑出来放到新队伍尾部

队伍一: 1, [3]

队伍二: 2, [6, 7]

新队伍: [1, 2]

第3步:队伍一的队头3较小,挑出来放到新队伍尾部

队伍一: 1, 3, []

队伍二: 2, [6, 7]

新队伍: [1, 2, 3]

第4步:队伍一已空,把队伍二剩下的部分按原顺序放到新队伍尾部

队伍一: 1, 3, []

队伍二: 2, 6, 7, []

新队伍: [1, 2, 3, 6, 7]

按照这个方法,我们把两个数组读进来之后,可以归并到一个新的数组里去,然后在新数组里计算中位数即可。

#include <cstdio>

int a[1000010], b[1000010], c[2000010];

int main()
{
	int m, n, half;
	scanf("%d %d", &m, &n);
	// 计算出元素总数 m+n 的一半,由于C++整数除法向下取整,以及数组从0开始计数的特点
	// 当 m+n 为奇数时,half 恰好值向中间元素的位置,例如 (2+3)/2=2,c[2]恰好是第3个元素
	// 当 m+n 为偶数时,half-1 和 half 恰为最中间的两个位置,例如 (2+2)/2=2,c[1]和c[2]位于正中位置
	half = (m + n) / 2;
	// 读入两个原数组 a 和 b
	for (int i = 0; i < m; ++i) scanf("%d", &a[i]);
	for (int i = 0; i < n; ++i) scanf("%d", &b[i]);
	// pa, pb: a 和 b 两个数组中还没有被归并部分的首位置;p: c 数组的尾部位置
	int pa = 0, pb = 0, p = 0;
	// 每次从 a 和 b 未归并部分的头部挑较小的那个数放到 c 的尾部,直到 a 或 b 中某一个数组被归并完
	while (pa < m && pb < n) a[pa] <= b[pb] ? c[p++] = a[pa++] : c[p++] = b[pb++];
	// 把剩余部分照搬到 c 的尾部完成归并
	while (pa < m) c[p++] = a[pa++];	// 如果 a 已经归并完了,这个循环不会进入
	while (pb < n) c[p++] = b[pb++];	// 同上,故两个循环只有一个会真正运行
	// 输出中位数
	printf("%.2lf\n", p % 2 ? (double) c[half] : (c[half] + c[half-1]) / 2.0);

	return 0;
}

这个算法思路很简单,程序很简洁,运行速度其实也挺快的,通常情况下这种方法就很好,推荐使用。但是归并算法有一个很大的问题是大量占用内存,因为它要把所有原数据都复制一份,所以空间复杂度高达 \(O(m+n)\)

归并的时间复杂度很显然是 \(O(m+n)\),所以这个基于归并的算法还远达不到不超过 \(O\bigl(\log(m+n)\bigr)\) 的时间复杂度要求。

我们需要一个更加快的算法,因为问题要求的时间复杂度为对数型,看到对数型的时间要求一定会想到二分法。这个世界上绝大多数达到对数时间的算法都是二分法。下面我们就要重点讲一讲怎么用二分法求两个有序数组的中位数。

3.4.10.2. 二分法求两个有序数组的中位数

根据中位数的定义,假设把 \(A,B\) 两个数组有序归并成数组 \(C\),取 \(half=\left\lfloor\frac{m + n}{2}\right\rfloor\),将 \(C\) 分为前后两半 \(C_l\)\(C_r\)。前半部分 \(C_l\)\(half\) 个数,后半部分 \(C_r\)\(m+n-half\) 个数。那么实际上我们已经找到中位数了:

如果 \(m+n\) 为奇数,中位数等于后半部分最小值 \(C_r[0]\)

如果 \(m+n\) 为偶数,则中位数为前半部分最大值和后半部分最小值的平均值 \(\frac{1}{2}\bigl(C_l[half-1]+C_r[0]\bigr)\)

但是现在我们不能真的去归并出数组 \(C\) 来,怎么办?我们当然还是要设法把所有的数划分成两半,并且在没有归并、没有排序的情况下判断是不是前一半中的最大值恰好小于等于后一半中的最小值,即 \(\max\{C_l\}\le \min\{C_r\}\)。请注意,现在 \(C_l\) 中的最大值不是 \(C_l[half-1]\) 了,\(C_r\) 中的最小值也不是 \(C_r[0]\) 了。如果找到了满足这个条件的一种划分,那么实际上中位数也就找到了:

如果 \(m+n\) 为奇数,中位数等于后半部分最小值 \(\min\{C_r\}\)

如果 \(m+n\) 为偶数,则中位数为前半部分最大值和后半部分最小值的平均值 \(\frac{1}{2}\bigl(\max\{C_l\}+\min\{C_r\}\bigr)\)

要找到这样一种划分,首先要有办法遍历所有可能的划分方法。

不失一般性,假设 \(m \le n\),即数组 \(A\) 的元素个数不多于数组 \(B\) 的元素个数。如果读入的数据情况相反,我们可以在编程的时候用一个小技巧把两个数组互换过来,这个技巧到后面再讲,现在先在这个假设下进行分析。

如前所述,\(C_l\) 中一定有一些数来自 \(A\),另一些来自 \(B\)\(C_r\) 也是如此。\(A\)\(B\) 都为 \(C\) 的前后两半贡献了一些数。不难想象,\(A\) 贡献给 \(C_l\) 的数一定都不大于它贡献给 \(C_r\) 的数,即 \(A\)\(C\) 前后两半的贡献其实也是它自己的一个前后划分。数组 \(B\) 当然也一样。所以,是 \(A\)\(B\) 各自的一个前后划分,组合成了 \(C\) 的前后两半。

于是我们可以先任意划分数组 \(A\),将 \(A[0:k]\) 贡献给 \(C_l\)\(A[k:m]\) 贡献给 \(C_r\),其中 \(0\le k\le m\)。当 \(k=0\) 时,\(A\) 中所有元素都贡献给了 \(C_r\);当 \(k=m\) 时则全部贡献给了 \(C_l\)

\(A\) 划分完毕后,实际上 \(B\) 的划分也就依此而确定了。这是因为 \(C_l\) 规定了是要有 \(half\) 个元素,所以既然 \(A\) 已经提供了 \(k\) 个,那么 \(B\) 就一定是提供 \(half-k\) 个。所以 \(B[0:half-k]\) 贡献给 \(C_l\),剩余的 \(B[half-k:n]\) 全部贡献给 \(C_r\)。因为我们已经限定了 \(m \le n\),所以确保了 \(B\) 的划分一定不会出界,即一定能确保 \(0\le half-k\le n\)。这是一个简单的不等式推导,请自行验证。

数组 \(A\) 一共有 \(m+1\) 种划分,其中必有一种是满足中位数划分的条件的,即 \(\max\{C_l\}\le \min\{C_r\}\),我们的任务于是就变成了找出这种划分。如果蛮力枚举,也已经把时间复杂度降到了 \(O\bigl(\min\{m,n\}\bigr)\)。如果对于任意一种划分 \(k\),我们有办法确定它是大了还是小了,即数组 \(A\)\(C_l\) 贡献其前 \(k\) 个元素是给多了还是给少了,那么我们就可以愉快地用二分法来查找正确的划分数 \(k\) 了,时间复杂度 \(O\bigl(\log\min\{m,n\}\bigr)\),问题完美解决。

要判断一个划分 \(k\) 是大了还是小了还是正好,一共有三种情况:

  1. \(1\le k \le m-1\),这时候 \(A[0:k]\) 贡献给了 \(C_l\)\(A[k:m]\)\(C_r\) 中。相应的,数组 \(B\) 的划分点为 \(j=half-k\),满足 \(1\le j\le n-1\),所以 \(B[0:j]\) 贡献给了 \(C_l\)\(C_r\) 中则是 \(B[j:n]\)。两个数组对前后两半都有贡献,如下所示:

                    C的前半部分     |     C的后半部分
    
    来自数组A:  A[0], ..., A[k-1]  |  A[k], ..., A[m-1]
    
    来自数组B:  B[0], ..., B[j-1]  |  B[j], ..., B[n-1]
    

    我们要比较前半部分的最大值和后半部分的最小值,前半部分的最大值是 \(A[k-1]\)\(B[j-1]\) 二者中的大者,后半部分的最小值是 \(A[k]\)\(B[j]\) 二者的小者。四者要相互比较大小,但是由于两个原数组本身是有序的,所以 \(A[k-1]\le A[k]\)\(B[j-1]\le B[j]\) 天生成立,它们不需要比较,我们只需要交叉对比即可:

    1. \(A[k-1]\gt B[j]\) 时,说明数组 \(A\) 贡献给 \(C_l\) 的数太多了,有过大的数混进前半部分去了,所以这时候我们应该减小 \(k\) 的值;

    2. \(B[j-1]\gt A[k]\) 时,说明数组 \(A\) 贡献给 \(C_l\) 的数太少了,有过小的数混进后半部分去了,所以这时候我们应该增大 \(k\) 的值;

    3. 其余情况,说明 \(A[k-1]\le B[j]\)\(B[j-1]\le A[k]\),这就是说 \(A[k-1]\)\(B[j-1]\) 二者中的大者小于等于 \(A[k]\)\(B[j]\) 二者的小者,即 \(C\) 前半部分的数都不大于后半部分的数,Bingo!我们找到了中位数。

  2. \(k=0\) 时,数组 \(A\) 中所有数都分在了 \(C_r\)\(C_l\) 全部由 \(B[0:j]\) 构成。此时 \(j=half\),当 \(m=n\)\(j=n\),会出现整个数组 \(B\) 都在 \(C_l\) 前半部分的情况,但是这没有关系,如下所示:

                    C的前半部分     |     C的后半部分
    
    来自数组A:                     |  A[0], ..., A[m-1]
    
    来自数组B:  B[0], ..., B[j-1]  |  B[j], ..., B[n-1]  (当 m == n 时,j == n,后半部分无 B 的元素)
    

    现在,数组 \(A\)\(C_l\) 没有贡献任何数,所以前半部分的最大值就是 \(B[j-1]\),而且 \(A\) 的划分数 \(k\) 也已经不能再减小了。因此只能有两种情况,要么 \(k\) 小了,要么这个划分正好。而 \(B\) 有没有对 \(C\) 的后半部分贡献元素是无关紧要的,因为要减小 \(k\) 的情况不存在,而且反正也没有数需要和它比大小:

    1. \(B[j-1]\gt A[k]\) 时,说明数组 \(A\) 贡献给 \(C_l\) 的数太少了,有过小的数混进后半部分去了,所以这时候我们应该增大 \(k\) 的值;

    2. 其余情况,说明这就是中位数的正确划分。

  3. \(k=m\) 时,数组 \(A\) 中所有数都被分在了 \(C_l\)\(C_r\) 全部由 \(B[j:n]\) 构成。此时 \(j=half-m\),当 \(m=n\)\(j=0\),会出现整个数组 \(B\) 都在 \(C_r\) 的情况,但是和情况2类似,这也没有关系,如下所示:

                    C的前半部分     |     C的后半部分
    
    来自数组A:  A[0], ..., A[m-1]  |
    
    来自数组B:  B[0], ..., B[j-1]  |  B[j], ..., B[n-1]  (当 m == n 时,j == 0,前半部分无 B 的元素)
    

    现在,数组 \(A\)\(C_r\) 没有贡献任何数,所以后半部分的最小值就是 \(B[j]\),而且 \(A\) 的划分数 \(k\) 也已经不能再增大了。因此只能有两种情况,要么 \(k\) 大了,要么这个划分正好。而 \(B\) 有没有对 \(C\) 的前半部分贡献元素是无关紧要的,因为要增大 \(k\) 的情况不存在,而且反正也没有数需要和它比大小:

    1. \(A[k-1]\gt B[j]\) 时,说明数组 \(A\) 贡献给 \(C_l\) 的数太多了,有过大的数混进前半部分去了,所以这时候我们应该减小 \(k\) 的值;

    2. 其余情况,说明这就是中位数的正确划分。

上面一共分析了三种情况,每种情况里又有2到3种小情况,显得比较复杂。但是我们发现所有七种小情况其实一共只有减小 \(k\)、增大 \(k\)、解是 \(k\) 这么三类,所以我们可以综合以上三种情况,整合成比较清晰的符合二分查找特征的三个分支:

  1. \(1\le k\le m\) 时,即当 \(k\gt0\) 时,如果 \(A[k-1]\gt B[j]\),那么减小 \(k\)

  2. \(0\le k\le m-1\) 时,即当 \(k\lt m\) 时,如果 \(B[j-1]\gt A[k]\),那么增大 \(k\)

  3. 其余情况,查找成功。

上面就是整个二分查找数组 \(A\) 的合适划分数 \(k\) 的过程。

这个二分算法,比我们前面学过的几个都要复杂,而且有以下几个不同之处:

  1. 问题必有解,所以二分查找的循环甚至可以用死循环,然后在判断找对答案的时候 break 循环。

  2. 找对答案的情况比没找对的情况要复杂,所以我们在编程的时候应该先对没找对答案的情况进行判断,最后剩余的就是找对了答案的情况。

  3. 由于必有解,所以在查找的边界处可能会有特殊情况需要照顾,在分析时要引起注意。

最后介绍一下怎样确保两个原数组一定是 \(m\le n\)。当然不是读完之后判断二者大小然后交换数组元素,那太笨了。我们可以定义一个二维数组用来保存两个原数组:

int d[2][1000010];

然后照常读入两个原数组的长度 mn,将数组 \(A\) 读入到 d[0] 中,将数组 \(B\) 读入到 d[1] 中。再定义两个整型变量:

int a = 0, b = 1;

接下来就是见证奇迹的时刻。我们比较两个原数组的长度,如果 m <= n 就什么都不做,如果 m > n 就交换 mn 的值,并且让 a = 1, b = 0。于是 d[a] 就是比较短的那个数组,长度为 md[b] 就是比较短的那个数组,长度为 n。搞定!这叫做“名实分离”,我说你是 \(A\) 你就是 \(A\),说你是 \(B\) 你就是 \(B\),不管你原本是谁。这个技巧在很多场景下有用,很实用。

下面是完整的程序:

#include <cstdio>

int d[2][1000010];
int a = 0, b = 1;

int main()
{
	int m, n, half;
	scanf("%d %d", &m, &n);
	for (int i = 0; i < m; ++i) scanf("%d", &d[a][i]);
	for (int i = 0; i < n; ++i) scanf("%d", &d[b][i]);
	half = (m + n) / 2;
	// 确保数组 d[a] 是比较短的那个数组,其长度为m,因此确保 m <= n
	if (m > n) { a = m; m = n; n = a; a = 1; b = 0; }
	
	// 二分查找数组 d[a] 分在前半部分的数字个数 k,查找范围为从 0 到 m,
	// 当 d[a] 中有 k 个数分在前半部分的时候,那么 d[b] 中必须有 j = half-k
	// 个数分在前半部分,这样前半部分一共 half 个数,后半部分 m+n-half 个。
	// 当 m+n 为偶数时,前后个数相等;奇数时,前半部分少一个。
	// 
	// 设前半部分最大值为 max_h,后半部分最小值为 min_r,那么当 m+n 为偶
	// 数时,中位数为 (max_h+min_r)/2.0;为奇数时,中位数为 min_r。
	//
	// 因为 m <= n,所以 m <= half <= n,当:
	//   k = 0 时:d[a] 的所有 m 个数都在后半部分
	//   k = m 时:d[a] 的所有 m 个数都在前半部分
	//   1 <= k <= m-1 时:d[a] 和 d[b] 各有一部分数在前后两个部分
	//
	// 又,各重要变量的取值范围如下:
	// 0 <= k <= m, m <= half <= n => 0 <= m-k <= j = half-k <= n-k <= n
	//
	// 此问题必有且仅有一个解,解的条件是:max_h <= min_r
	//
	// 查找时应判断 d[a] 提供给前半部分的数字个数 k 是多了还是少了:
	//     如果 d[a] 贡献给前半部分的最后一个数字 d[a][k-1] 比 d[b] 留在
	//         后半部分的第一个数字 d[b][j] 大,那么说明 k 多了
	//     如果 d[b] 贡献给前半部分的最后一个数字 d[b][j-1] 比 d[a]
	//         留在后半部分的第一个数字 d[a][k] 大,那么说明 k 少了
	// 所以,进行下面的分析:
	//     当 1 <= k <= m-1 时:
	//        IF d[a][k-1] > d[b][j] THEN right = k - 1
	//        IF d[b][j-1] > d[a][k] THEN left = k + 1
	//     当 k == 0 时:
	//        这时已经不可能再减少 k 了,也不存在 d[a][k-1],只需要判断
	//        IF d[b][j-1] > d[a][k] THEN left = k + 1
	//     当 k == m 时:
	//        这是已经不可能再增加 k 了,也不存在 d[a][k],只需要判断
	//        IF d[a][k-1] > d[b][j] THEN right = k - 1
	//     其他情况都说明查找成功。
	// 以上三种情况可以合并为两种:
	//     IF k < m && d[b][j-1] > d[a][k] THEN left = k + 1
	//     IF k > 0 && d[a][k-1] > d[b][j] THEN right = k - 1
	// 剩余情况说明查找成功,计算 max_h 和 min_r 并得出中位数即可
	//
	int left = 0, right = m, k, j;
	// d[a], d[b] 各自在前半部分的最大值和后半部分的最小值
	int max_h, min_r;
	while (left <= right) {
		k = (left + right) / 2;
		j = half - k;
		if (k < m && d[b][j-1] > d[a][k])
			left = k + 1;	// k 小了
		else if (k > 0 && d[a][k-1] > d[b][j])
			right = k - 1;	// k 大了
		else {			// Bingo! 找到中位数的正确划分
			// 计算左半部分的最大值 max_h
			if (k == 0)
				max_h = d[b][j-1];
			else if (j == 0)
				max_h = d[a][k-1];
			else
				max_h = d[a][k-1] > d[b][j-1] ? d[a][k-1] : d[b][j-1];
			// 计算右半部分的最小值 min_r
			if (k == m)
				min_r = d[b][j];
			else if (j == n)
				min_r = d[a][k];
			else
				min_r = d[a][k] < d[b][j] ? d[a][k] : d[b][j];
			// 退出循环
			break;
		}
	}
	printf("%.2lf\n", (m + n) % 2 ? (double) min_r : (max_h + min_r) / 2.0);

	return 0;
}

最后留一个思考问题:为什么不需要在开始查找之前对“有一个数组为空”的情况进行特判?