欧式距离如何应对缺失值
biobabble作者群有一个小伙伴提出来的问题。
没错，凡是给biobabble投稿的小伙伴都会被拉入作者群。
话说有一个矩阵，我们让它某个值缺失，算出来的距离和原来的竟然是一样的：
> dat = matrix(1:10, 2)
> dist(dat)
1
2 2.236068
> dat[1,1] = NA
> dist(dat)
1
2 2.236068
这有点反直觉，实际上你再搞个缺失值，它算出来还是一样的：
> dat[1,3] = NA
> dist(dat)
1
2 2.236068
但以我的经验，我通常不敢轻易说人家有bug，我们得确认一下。
于是我读了R的源代码中欧式距离的代码：
static double R_euclidean(double *x, int nr, int nc, int i1, int i2)
{
double dev, dist;
int count, j;
count= 0;
dist = 0;
for(j = 0 ; j < nc ; j++) {
if(both_non_NA(x[i1], x[i2])) {
dev = (x[i1] - x[i2]);
if(!ISNAN(dev)) {
dist += dev * dev;
count++;
}
}
i1 += nr;
i2 += nr;
}
if(count == 0) return NA_REAL;
if(count != nc) dist /= ((double)count/nc);
return sqrt(dist);
}
把NA干掉？
首先从代码上看，缺失值会被去掉，这就是count计数看有多少有效值的作用。如果我们简单粗爆地把NA值去掉，我们会算出来什么值呢？
- 完整数据的情况下，距离是2.236068
- dat[1,1] = NA的情况下，距离是2
- dat[1,c(1,3)] = NA的情况下，距离是1.732051
从这里我们可以看到，缺失值越多，距离是越小的，这个从欧式距离的公式上也反映出来，所有的量都是正值，少了就小了。所以如果我们简单地去NA值，那么对距离的估算会低估。
怎么补救？
这就是倒数第二句干的事情：
if(count != nc) dist /= ((double)count/nc);
如果有缺少，用已知的点的均值来估计。也就是用mean((qi-pi)^2)去代替缺失点。
不难验证一下，前面第二种情况时距离是2，那么被一下缺失值：
> sqrt(2^2 / (4/5))
[1] 2.236068
出来的就是第一种情况时的值了。
如果你看不明白的话，我再拆一下给你看：
- 2^2 是总和
- 2^2⁄4 是均值
- 2^2 + 2^2⁄4 就是补了缺失值的总和
- 2^2 / (4⁄5)等同于上面的值
- sqrt(2^2 / (4⁄5))算出来就是距离，请对照公式来看一下。
为什么能够完美估计原来的值？
这也是一开头让大家疑惑的地方，为什么做为估计出来的值能和真实的值一样？
这其实是个巧合，真的是巧合。你如果是一般的随机数拿来试试，肯定无法一模一样的值。
当然对于一个determinant的计算，巧合也是算出来的，必定是巧合，原因必须存在于数值之中。
如果你对dat算相关系数的话，就明白了。我们对原始值算一下：
> cor(dat[1,], dat[2,])
[1] 1
能够完美估算出距离，其实等同于能够完美复原缺失点，而这种情况只在于数据的相关系数等于一的情况下，数据的相关性越好，也就能够复原得越好，这很好理解。
欧式距离对于outlier是很敏感的，你可以试试不同的值去替换掉原始值，再看看缺失值的情况下，计算出来的距离和不缺失的情况下，差距有多少。
写到这里，我想到了欧式距离可以拿来检测outlier啊。比如我们每次让一个点缺失，比如说我们有100个点，就可以计算出100个距离，这些距离各不相等（正常情况下数据相关系统不会刚刚好是1），它们和真实值的差值也可以算出来，这100个差值里面，差值较大的，相对应的缺失点，就是outlier。
