Categories
不学无术

LeetCode 287. Find the Duplicate Number | 智商被碾压!

看了LeetCode上牛人的解答,感觉我的智商被碾压了~

原题

Given an array nums containing n + 1 integers where each integer is between 1 and n (inclusive), prove that at least one duplicate number must exist. Assume that there is only one duplicate number, find the duplicate one.
给定一个包含n+1个整数的数组nums,该数组中的每个整数都在[1,n]的范围内,假设数组中只有一个数是重复的,找到那个重复的数。
Note|提示:

  1. You must not modify the array (assume the array is read only).|    不能改动数组内容
  2. You must use only constant, O(1) extra space.|    在O(1)的空间复杂度内完成
  3. Your runtime complexity should be less than O(n2).|    时间复杂度不超过O(n^2)
  4. There is only one duplicate number in the array, but it could be repeated more than once.|    只有一个重复数字,但是可能重复超过一次

思路1:O(nlogn)时间复杂度

英文版解释:https://leetcode.com/discuss/60830/solutions-explanation-space-without-changing-input-array
这个解法涉及到鸽巢原理。
假设有n=10的一个数组(大小是n+1),那我目前的搜索范围是[1,10],先找中间的数mid=5.现在我们遍历数组的所有元素并统计<=5的元素个数,记作n_lteq5好了,那么有:

  • 如果n_lteq5>5,那就有6个数字占了本来只有5个的坑位,那目标数字肯定是在[1,5]的范围内了;
  • 如果n_lteq5<=5,那前面的元素都挺守规矩的,得看看[6,10]里头哪个数字在作怪;

这样每一次判断,问题的规模都会缩小一半,一共n+1个数字,时间复杂度是O(nlogn)
感觉碰到了很多二分思想的题目啊,关键点就是要找到一个什么划分方式把问题的规模缩小掉,而在这道题里头用到的就是鸽巢原理。

思路1代码

class Solution {
public:
    int countLtEq(const vector<int>& nums, int n){
        int count = 0;
        for(int i=0; i<nums.size(); i++){
            if(nums[i] <= n)
                count ++;
        }
        return count;
    }
    int findDuplicate(vector<int>& nums) {
        int n = nums.size() - 1;
        int min = 1, max = n;
        int mid;
        int ltEqCnt;
        while(min < max){
            mid = (min + max)/2;
            ltEqCnt = countLtEq(nums, mid);
            if(ltEqCnt <= mid){
                //[min, mid]是正常的,找另一边
                min = mid+1;
            } else {
                max = mid;
            }
        }
        return min; // min == max
    }
};

 

插曲:如果可以改变数组元素呢?

感觉在剑指Offer上遇到过类似的题目。如果可以改变数组元素的话,可以按照下面的思路进行,时间复杂度是O(n)。基本思想是,一个萝卜一个坑,把搜索到的数放到自己应该在的位置上,比如1应该放在第一个位置上,如果后面再找到一个1,我们就会发现1这个位置的坑被前一个1占了,那就冲突了。
方便起见,数组的位置按1,2,3,…描述(而非0,1,2…)。假设有数组[5,4,2,3,1,5,6] ,当前工作的位置pos=1

  • 先看第pos=1个数字5,在设想的排好序的序列里,5应该就在第5个位置上,因此我们把5和第5个位置上现有的数字1进行调换,结果就变成了[1,4,2,3,5,5,6] ;
  • 看换回来的数字1,由于已经在1号位上,这次处理就结束了,往前挪一下,pos=2,看下一个数字;
  • pos=2的数字是4,本来应该在4号位上,因此把4和4号位上的 3调换,数组变成[1,3,2,4,5,5,6] ;
  • 这时候事情没完,因为换回来的3不在应该有的3号位上,因此把3和3号位上占着的2调换,数组变成[1,2,3,4,5,5,6] ;
  • 换回的数字2是正确的位置,因此pos++,变成pos=3,类似地检查pos=4, 5发现都符合要求;
  • pos=6时,此时位置上的元素是5,于是我们想把这个5归为到5号位,但是检查5号位的时候我们发现,上面已经有一个5了,这就说明5是重复数值了,搞定

可能有问,那如果原来的数组是[5,4,2,3,1,6,5] 呢?很显然,如果前面找到pos=6的时候都没发生冲突,那罪魁祸首肯定在最后个元素上面了。
这种方法只要遍历一次数组并做相应替换就行,除了会改变数组内容外,其他条件都满足,可惜啊。

思路2:O(n)时间复杂度

大神,大家快去膜拜:http://keithschwarz.com/interesting/code/?dir=find-duplicate
据说越是经典的代码行数越是写的简单,然而你并不能看懂…
似乎是把数组中的值看做单链表的next指针的值的,还在理解中~这个算法理论上是正确的。
 

Categories
不学无术

(转载)从最大似然到EM算法浅解

感觉是我看到解释得最容易理解的一个版本
http://blog.csdn.net/zouxy09/article/details/8537620

Categories
不学无术

LeetCode #4. Median of Two Sorted Arrays

题目链接https://leetcode.com/problems/median-of-two-sorted-arrays/
参考文献

思路:
实在是脑子转不过来,看了一下解法,然后就把自己的理解讲一下吧。
有两个规律需要说下:
1.在某个数列头尾删去相同数量的元素,其中位数是不变的。比如[1,2,3,4]删去1和4.
2.如果两个数列的中位数是m1, m2,那么合并后的数列的中位数的值肯定在[m1,m2]范围内。
既然题目是要Log级别的复杂度,就要想到要用折半法把问题的规模降低…举个例子,假设有[1,2,4,6,6,8]和[1,5,7,8,10,25,36]两个数列,很容易算出他们的(前序 )中位数分别是4和8,而他们合并后的中位数的值在[4,8]之间。
然后可以删去无用的值,[1,2]与[10,25,36]. 但由于规律1的限制,如果后面删了3个的话,合并后删减的数列会发生变化(本例中奇偶性都变了),我们只能删去min(2,3)=2个元素…这种删值从理论上讲就是折半了,然后搞一波递归就行了。
最最最基本的思路如下:

两个有序数列Ary1, Ary2(假设升序排列)的中位数分别为m1,m2
if m1 == m2:
    LUCKY!!!
elif m1 < m2:
    取(Ary1的后半段,Ary2的前半段)再算
else
    取(Ary1的前半段,Ary2的后半段)再算

由于题设里面俩数组是不等长的,所以里对其中某个数组n<=2时处理方法需要注意…
代码:

int max(int x, int y) {
    return x>y ? x : y;
}
int min(int x, int y) {
    return x<y ? x : y;
}
double getMedian(int* nums, int size) {
	if (size == 0)
		return -1;
	if (size % 2 == 0) {
		//Even number
		return (nums[size / 2] + nums[size / 2 - 1]) / 2.0;
	}
	else {
		return nums[(size - 1) / 2];
	}
}
double findMedianMerge(int* nums1, int nums1Size, int* nums2, int nums2Size) {
	int count = 0;
	int n1Count = 0;
	int totalSize = nums1Size + nums2Size;
	int t1, t2;
	if (totalSize % 2 == 0) {
		t1 = totalSize / 2;
		t2 = t1 + 1;
	}
	else {
		t1 = -1;
		t2 = totalSize / 2 + 1;
	}
	double sum = 0;
	int curVal;
	while (n1Count < nums1Size) {
		if (*nums1 < *nums2) {
			n1Count++;
			curVal = *nums1;
			nums1++;
		}
		else {
			curVal = *nums2;
			nums2++;
		}
		count++;
		if (count == t1 || count == t2) {
			sum += curVal;
		}
		else if (count >= t2) {
			break;
		}
	}
	if (count < t1) {
		//Not finished
		nums2 += (t1 - count - 1);
		count = t1;
		sum += *nums2;
		nums2++;
	}
	if (count < t2) {
		nums2 += (t2 - count - 1);
		sum += *nums2;
	}
	if (t1 == -1)
		return sum;
	else
		return sum / 2.0;
}
double findMedianSortedArrays(int* nums1, int nums1Size, int* nums2, int nums2Size) {
	int totalSize = nums1Size + nums2Size;
	if (nums1Size > nums2Size) {
		//Keep ary1 shorter than ary2
		int* temp = nums1;
		int t = nums1Size;
		nums1 = nums2;
		nums1Size = nums2Size;
		nums2 = temp;
		nums2Size = t;
	}
	if (nums1Size == 0) {
		return getMedian(nums2, nums2Size);
	} else if (nums1Size == 1) {
		if (nums2Size == 1) {
			return (nums1[0] + nums2[0]) / 2.0;
		}
		else {
			return findMedianMerge(nums1, nums1Size, nums2, nums2Size);
		}
	}
	else if (nums1Size == 2) {
		if (nums2Size == 2) {
			return (max(nums1[0], nums2[0]) + min(nums1[1], nums2[1])) / 2.0;
		}
		else {
			return findMedianMerge(nums1, nums1Size, nums2, nums2Size);
		}
	}
	//Slice
	double m1 = getMedian(nums1, nums1Size);
	double m2 = getMedian(nums2, nums2Size);
	if (m1 == m2) {
		//Got lucky
		return m1;
	}
	else if (m1 < m2) {
		//Slice L-nums1 and R-nums2
		int cut = (nums1Size - 1) / 2;
		return findMedianSortedArrays(nums1 + cut, nums1Size - cut, nums2, nums2Size - cut);
	}
	else {
		//Slice R-nums1 and L-nums2
		int cut = (nums1Size - 1) / 2;
		return findMedianSortedArrays(nums1, nums1Size - cut, nums2 + cut, nums2Size - cut);
	}
}

 

Categories
不学无术

Glickman 的ELO算法

Glickman的ELO算法被搬到众包中实现,可以参考这篇

Bashir, M., Anderton, J., Wu, J., Golbus, P. B., Pavlu, V., & Aslam, J. A. (2013). A Document Rating System for Preference Judgements. In Proceedings of the 36th International ACM SIGIR Conference on Research and Development in Information Retrieval (pp. 909–912). New York, NY, USA: ACM. http://doi.org/10.1145/2484028.2484170
但是文章里算法描述的部分似乎有写错的地方,可能是作者不小心,把公式(7)的分母部分写错了(原文中的字母j没有换成自己的B)~
Glickman的原文是:
Glickman, M. E. (1999). Parameter estimation in large dynamic paired comparison experiments. Journal of the Royal Statistical Society: Series C (Applied Statistics), 48(3), 377-394.
QQ截图20150622165240
其中m为对手个数, 为与对手B进行比赛的次数; 为本局比赛A的结果(1-胜,0-负)
F是一个常量,(Bashir et al. 2013)中F=200, (Glickman, 1998)中F=400.
 
The rating algorithm is implemented as follows.
(a) Collect game outcome data over a rating period.
(b) At the end of the period, update players’ rating distributions due to game outcomes from their preperiod (prior) rating distributions.
(c) Subsequently update players’ rating distributions due to the passage of time.
Categories
不学无术

组合 递归

荒废计算机太久了,以至于这么简单一个问题想了我大半天,其实就是数组求组合输出的方法,写下此日志提醒自己不能忘了老本!
题目的话大概可以这样表述

给出一个数据A={a_0,a_1,a_2…a_n}(其中n可变),打印出该数值元素的所有组合

需要使用递归的方法解决比较好,分而治之

给一大团数据{a0,a1, a2, …, an}:
1. 如果就剩一个元素的话,我就输出它,结束;
2.不输出a0,把除了我之外的余下一团交由递归处理;
3.要输出a0自己,然后余下一团交由递归处理;

其中case 3里面需要注意下一个递归要记住a0是要输出的,所以实践起来就是带了个前缀一样的东西记住。感觉这种排列组合的问题应该也可以用栈来解决,但归根结底还是分治。
用Python随意写了一个:

# coding:utf-8
def foo(inary, appendix=''):
    if len(inary) == 0:
        return
    #Case 0
    print appendix + str(inary[0])
    if len(inary) < 2:
        return
    #Case 1
    foo(inary[1:], appendix)
    #Case 2
    #print str(inary[0]),
    foo(inary[1:], appendix+str(inary[0]))

测试输出:

>>> foo([1, 2, 3, 4])
1
2
3
4
34
23
24
234
12
13
14
134
123
124
1234

 

Categories
不学无术

KM算法 详解+模板

本文转载自:http://blog.sina.com.cn/s/blog_691ce2b701016reh.html
 
先说KM算法求二分图的最佳匹配思想,再详讲KM的实现。
【KM算法求二分图的最佳匹配思想】

对于具有二部划分( V1, V2 )的加权完全二分图,其中 V1= { x1, x2, x3, … , xn }, V2= { y1, y2, y3, … , yn },边< xi, yj >具有权值 Wi,j 。该带权二分图中一个总权值最大的完美匹配,称之为最佳匹配。
 
记 L(x) 表示结点 x 的标记量,如果对于二部图中的任何边<x,y>,都有 L(x)+ L(y)>= Wx,y,我们称 L 为二部图的可行顶标。
设 G(V,E) 为二部图, G'(V,E’) 为二部图的子图。如果对于 G’ 中的任何边<x,y> 满足, L(x)+ L(y)== Wx,y,我们称 G'(V,E’) 为 G(V,E) 的等价子图。
 
定理一:设 L 是二部图 G 的可行顶标。若 L 等价子图 G有完美匹配 M,则 M 是 G 的最佳匹配。
证明:由于 GL 是 G 的等价子图,M 是 GL 的完美匹配,所以,M 也是 G  的完美匹配。以由于对于匹配 M 的每条边 e ,都有 e∈ E( GL ),而且 M 中每条边覆盖每个顶点正好一次,所以
W( M )= å W(e), e∈ M = å L(x), x∈ V
另一方面,对于 G 的任何完美匹配 M’ 有
W( M’ )= å W(e), e∈ M’ <= å L(x), x∈ V
于是 W( M )>= W( M’ ),即 M 是 G 的最优匹配。
 
由上述定理,我们可以通过来不断修改可行顶标,得到等价子图,从而求出最佳匹配。
就像匈牙利算法一样,我们依次为每一个顶点 i 寻找增广路径,如果寻找增广路径失败,我们就修改相应的可行顶标,来得到增广路径。
如图:
 1  2  3  |
 3  2  4  |
 2  3  5  |
若要对这个完全二分图求最佳匹配
 
初始化:
Lx(1)= max{ y| w(1,y), 1<= y<= 3 }= max{ 1, 2, 3 }= 3, Ly(1)= 0
Lx(2)= max{ 3, 2, 4 }= 4, Ly(2)= 0
Lx(3)= max{ 2, 3, 5 }= 5, Ly(3)= 0;
我们建立等价子图( 满足 Lx(x)+ Ly(y)== W(x,y) ) 如下:
km算法求二分图最佳匹配
对于该图,运用匈牙利算法对 X 部顶点 1 求增广路径,得到一个匹配,如图( 红色代表匹配边 ):km算法求二分图最佳匹配
 对 X 部顶点 2 求增广路径失败,寻找增广路径的过程为 X 2-> Y 3-> X 1。我们把寻找增广路径失败的 DFS 的交错树中,在 X 部顶点集称之为 S, 在 Y 部的顶点集称之为 T。则 S= { 1, 2 },T= { 3 }。现在我们就通过修改顶标值来扩大等价子图,如何修改。
 
1)   我们寻找一个 d 值,使得 d= min{ (x,y)| Lx(x)+ Ly(y)- W(x,y), x∈ S, y∉ T },因些,这时 d= min{
Lx(1)+Ly(1)-W(1,1),  Lx(1)+Ly(2)-W(1,2),  Lx(2)+Ly(1)-W(2,1),  Lx(2)+Ly(2)-W(2,2) }=
min{ 3+0- 1, 3+0-2,  4+0-3,  4+0-2 }= min{ 2, 1, 1, 2 }= 1。
寻找最小的 d 是为了保证修改后仍满足性质对于边 <x,y> 有 Lx(x)+ Ly(y)>= W(x,y)。
 
2)   然后对于顶点 x
1. 如果 x∈ S 则 Lx(x)= Lx(x)- d。
2. 如果 x∈ T 则 Ly(x)= Ly(x)+ d。
3. 其它情况保持不变。
如此修改后,我们发现对于边<x,y>,顶标 Lx(x)+ Ly(y) 的值为
1.  Lx(x)- d+ Ly(y)+ d,  x∈ S, y∈ T。
2.  Lx(x)+ Ly(y),  x∉ S,  y∉ T。
3.  Lx(x)- d+ Ly(y), x∈ S, y∉ T。
4.  Lx(x)+ Ly(y)+ d, x∉ S,  y∈ T。
易知,修改后对于任何边仍满足 Lx(x)+ Ly(y)>= W(x,y),并且第三种情况顶标值减少了 d,如此定会使等价子图扩大。
 
就上例而言: 修改后 Lx(1)= 2, Lx(2)= 3, Lx(3)= 5, Ly(1)= 0, Ly(1)= 0, Ly(2)= 0, Ly(3)= 1。
这时 Lx(2)+Ly(1)=3+0=3= W(2,1),在等价子图中增加了一条边,等价子图变为:
 km算法求二分图最佳匹配
如此按以上方法,得到等价子图的完美匹配。
 
另外计算 d 值的时候可以进行一些优化。
定义 slack(y)= min{ (x,y)| Lx(x)+ Ly(y)- W(x,y),x∈ S,  y∉ T }
这样能在寻找增广路径的时候就顺便将 slack 求出。

(以上为摘上网络)
【KM算法及其具体过程】
(1)可行点标:每个点有一个标号,记lx[i]为X方点i的标号,ly[j]为Y方点j的标号。如果对于图中的任意边(i, j, W)都有lx[i]+ly[j]>=W,则这一组点标是可行的。特别地,对于lx[i]+ly[j]=W的边(i, j, W),称为可行边
(2)KM 算法的核心思想就是通过修改某些点的标号(但要满足点标始终是可行的),不断增加图中的可行边总数,直到图中存在仅由可行边组成的完全匹配为止,此时这个 匹配一定是最佳的(因为由可行点标的的定义,图中的任意一个完全匹配,其边权总和均不大于所有点的标号之和,而仅由可行边组成的完全匹配的边权总和等于所 有点的标号之和,故这个匹配是最佳的)。一开始,求出每个点的初始标号:lx[i]=max{e.W|e.x=i}(即每个X方点的初始标号为与这个X方 点相关联的权值最大的边的权值),ly[j]=0(即每个Y方点的初始标号为0)。这个初始点标显然是可行的,并且,与任意一个X方点关联的边中至少有一条可行边
(3)然后,从每个X方点开始DFS增广。DFS增广的过程与最大匹配的Hungary算法基本相同,只是要注意两点:一是只找可行边,二是要把搜索过程中遍历到的X方点全部记下来(可以用vst搞一下),以进行后面的修改;
(4) 增广的结果有两种:若成功(找到了增广轨),则该点增广完成,进入下一个点的增广。若失败(没有找到增广轨),则需要改变一些点的标号,使得图中可行边的 数量增加。方法为:将所有在增广轨中(就是在增广过程中遍历到)的X方点的标号全部减去一个常数d,所有在增广轨中的Y方点的标号全部加上一个常数d,则 对于图中的任意一条边(i, j, W)(i为X方点,j为Y方点):
<1>i和j都在增广轨中:此时边(i, j)的(lx[i]+ly[j])值不变,也就是这条边的可行性不变(原来是可行边则现在仍是,原来不是则现在仍不是);
<2>i在增广轨中而j不在:此时边(i, j)的(lx[i]+ly[j])的值减少了d,也就是原来这条边不是可行边(否则j就会被遍历到了),而现在可能是;
<3>j在增广轨中而i不在:此时边(i, j)的(lx[i]+ly[j])的值增加了d,也就是原来这条边不是可行边(若这条边是可行边,则在遍历到j时会紧接着执行DFS(i),此时i就会被遍历到),现在仍不是;
<4>i和j都不在增广轨中:此时边(i, j)的(lx[i]+ly[j])值不变,也就是这条边的可行性不变。
这 样,在进行了这一步修改操作后,图中原来的可行边仍可行,而原来不可行的边现在则可能变为可行边。那么d的值应取多少?显然,整个点标不能失去可行性,也 就是对于上述的第<2>类边,其lx[i]+ly[j]>=W这一性质不能被改变,故取所有第<2>类边的 (lx[i]+ly[j]-W)的最小值作为d值即可。这样一方面可以保证点标的可行性,另一方面,经过这一步后,图中至少会增加一条可行边。
(5)修改后,继续对这个X方点DFS增广,若还失败则继续修改,直到成功为止;
(6)以上就是KM算法的基本思路。但是朴素的实现方法,时间复杂度为O(n4)——需要找O(n)次增广路,每次增广最多需要修改O(n)次顶标,每次修改顶 标时由于要枚举边来求d值,复杂度为O(n2)。实际上KM算法的复杂度是可以做到O(n3)的。我们给每个Y顶点一个“松弛量”函数slack,每次开 始找增广路时初始化为无穷大。在寻找增广路的过程中,检查边(i,j)时,如果它不在相等子图中,则让slack[j]变成原值与 A[i]+B[j]-w[i,j]的较小值。这样,在修改顶标时,取所有不在交错树中的Y顶点的slack值中的最小值作为d值即可。但还要注意一点:修 改顶标后,要把所有不在交错树中的Y顶点的slack值都减去d。
【求二分图的最小匹配】
只需把权值取反,变为负的,再用KM算出最大权匹配,取反则为其最小权匹配。

hdoj 2255
#include 
#include 
#define M 310
#define inf 0x3f3f3f3f
int n,nx,ny;
int link[M],lx[M],ly[M],slack[M];    //lx,ly为顶标,nx,ny分别为x点集y点集的个数
int visx[M],visy[M],w[M][M];
int DFS(int x)
{
    visx[x] = 1;
    for (int y = 1;y <= ny;y ++)
    {
        if (visy[y])
            continue;
        int t = lx[x] + ly[y] - w[x][y];
        if (t == 0)       //
        {
            visy[y] = 1;
            if (link[y] == -1||DFS(link[y]))
            {
                link[y] = x;
                return 1;
            }
        }
        else if (slack[y] > t)  //不在相等子图中slack 取最小的
            slack[y] = t;
    }
    return 0;
}
int KM()
{
    int i,j;
    memset (link,-1,sizeof(link));
    memset (ly,0,sizeof(ly));
    for (i = 1;i <= nx;i ++)            //lx初始化为与它关联边中最大的
        for (j = 1,lx[i] = -inf;j <= ny;j ++)
            if (w[i][j] > lx[i])
                lx[i] = w[i][j];
    for (int x = 1;x <= nx;x ++)
    {
        for (i = 1;i <= ny;i ++)
            slack[i] = inf;
        while (1)
        {
            memset (visx,0,sizeof(visx));
            memset (visy,0,sizeof(visy));
            if (DFS(x))     //若成功(找到了增广轨),则该点增广完成,进入下一个点的增广
                break;  //若失败(没有找到增广轨),则需要改变一些点的标号,使得图中可行边的数量增加。
                        //方法为:将所有在增广轨中(就是在增广过程中遍历到)的X方点的标号全部减去一个常数d,
                        //所有在增广轨中的Y方点的标号全部加上一个常数d
            int d = inf;
            for (i = 1;i <= ny;i ++)
                if (!visy[i]&&d > slack[i])
                    d = slack[i];
            for (i = 1;i <= nx;i ++)
                if (visx[i])
                    lx[i] -= d;
            for (i = 1;i <= ny;i ++)  //修改顶标后,要把所有不在交错树中的Y顶点的slack值都减去d
                if (visy[i])
                    ly[i] += d;
                else
                    slack[i] -= d;
        }
    }
    int res = 0;
    for (i = 1;i <= ny;i ++)
        if (link[i] > -1)
            res += w[link[i]][i];
    return res;
}
int main ()
{
    int i,j;
    while (scanf ("%d",&n)!=EOF)
    {
        nx = ny = n;
      //  memset (w,0,sizeof(w));
        for (i = 1;i <= n;i ++)
            for (j = 1;j <= n;j ++)
                scanf ("%d",&w[i][j]);
        int ans = KM();
        printf ("%dn",ans);
    }
    return 0;
}
Categories
不学无术 木有技术

协同过滤文章两篇

探索推荐引擎内部的秘密,第 2 部分_ 深入推荐引擎相关算法 – 协同过滤 协同过滤

Categories
不学无术

interval scale和ratio scale的区别

interval scale 等距尺度:又称区间尺度,衡量的数字可代表大小或优劣顺序,数字的差距直接表示其差异的程度。
ratio scale 比例尺度:物与物的相对比较,表明各种物体间的相对度量关系。
参考百度百科

Categories
不学无术 木有技术

浅谈PCA 人脸识别

版权声明:转载时请以超链接形式标明文章原始出处和作者信息及本声明
http://leen2010.blogbus.com/logs/124631640.html

前几天讨论班我讲了基于PCA的人脸识别,当时我自己其实也只是知道这个算法流程,然后基于该算法利用c++实现了,效果还不错。后来跟师兄一起讨论的时候,才发现这个PCA还是有相当深刻的意义。
PCA的算法:矩阵C=AAT,A的每一列是一张人脸(将一张人脸图片用一个列向量表示,即对于128*128的图片,将视为16384维的列向量),A的列数就是图片的张数。算法就是求矩阵C的特征向量,每个向量称之为特征脸[1]。为了简单,只取其中部分的特征向量,这些特征向量对应于某些特征值,通常是前M个大的特征值。这样便得到了M个特征向量。接下来就是将每张图片在这M个特征向量上做投影,得到一个M维的权重向量w=(w1,…wM),一个人可能有多张图片,于是将对应于这个人的权重向量做一个平均作为这个人的权重向量。然后对于每个新来的人脸,先求得一个权重向量,然后与人脸库中每个人的权重向量做比较,如果小于某个阈值,则认为他是人脸库中的这个人;否则视为unknown。当然,文章中还给出了另外一个判断一张图像是否是人脸的方法,这里不再讨论。对于计算的时候我们实际上求的是ATA,至于二者的关系,可以参考文章[1],因为与后面讨论的关系不大,这里不细说了。
有个上述简介,我们就大概知道了PCA到底是怎么做的了。刨根问底一下,C到底是什么,C的特征向量又到底是什么,对应于特征值大的特征向量有有什么意义。
1.C到底是什么?我们先看一下C的(i,j)元素是什么。很简单,就是A的第i行与AT的第j列的乘积。那么A的第i行又是什么呢?就是每张图片的第i个像素构成的一个向量。则C的(i,j)元素就是每张图片的第i个像素构成的向量与每张图片的第j个像素构成的向量的乘积。回忆一下概率论中的协方差cov(x,y)=E((x-x)(y- y)),这里x– 是x的平均。如果我们把图片的第i个像素视为一个随机变量Xi,那么人脸库中的每张人脸的第i个像素的取值即为这个随机变量的所有取值。根据注,A的第i行的每个值都已经减去了Xi的均值,因此C的(i,j)元素就是随机变量Xi与Xj的协方差。到此,我们已经知道C是什么了。
2.C的特征向量是什么。我们知道对C求特征向量是找到一个可逆矩阵Q,使得Q-1CQ是一个对角阵D,D的元素即为特征值,Q中的每一列即为特征向量,与特征值所在的列一一对应。注意,因为C是实对称阵,故必然可以对角化。由于C是对称的,C的特征向量是正交的,因此Q便是一个正交阵,故Q-1即为QT。先从简单的角度来看,假设C已经是一个对角阵了,并且对角元素依次递减。即随机变量Xi与Xj(i!=j)时是不相关的,而Xi的方差则随着i的增大而减小。也就是前几个像素是方差比较大的像素,即在第一个像素上,每张图片的值变化最大,第二个像素次之,依此类推。举个例子,假设第一个像素表示人脸的额头上的某个点(也许你会问,第一个像素不是图片的最左上角的那个像素吗?为什么会是额头上的某个点,后面会说明为什么可以这么假设),而在这个点上,有些人有痣,有些人没有,那么这一点的方差就会比相邻的额头上的其他点大,因为其他点上,每个人的额头都差不多,因此其像素值也差不多,方差就小了。现在来考虑QTAATQ=D,这里依然假设D中的元素依次递减。对于QTA,QT的每一行是一个特征向量,QT的第i行与A的每一列相乘得到QTA的每一行,从矩阵的角度也可以看作是用QT对A做一个变换。记住这里的QTA可以看做前面讨论的A,那么变换的结果就是使得QTA前面的行对应的是方差大的像素,而且这个变换会把方差最大的放到了第一行(因为D的降序排列),这里也就解释了为什么前面那个例子可以认为第一个像素是额头上的某个点,因为做了变换。我们选择了QT的前M行作为特征向量,因为他们对应了前M个大的特征值。这里可以举一个直观但是不一定恰当的例子,人的头发部分对应的那些像素,经过QT变换后回到某个像素上,那么这个像素会是QTA中的哪个位置呢,我认为应该是QTA的列中靠下面的像素,因为在头发这个像素的地方每个人基本都是头发(这句话很别扭,我是想说,对于比较正规的人脸库,即每张人脸不会有太大的变化,某个人头发的对应的那几个像素对于其他人来说也都是头发,因此变换后这个像素对于每个人都差不多,方差小,固然会在比较靠后的位置了)。QTA的每一列的前M个元素对应的就是每张人脸的权重向量w。因此每张人脸的权重向量的同一个分量可以看作是新的像素,这些新的像素对应着方差依次减小的像素。对于一张新来的人脸,让他在特征向量上作投影得到权重向量,在这些方差大的像素处,如果跟某个人的比较接近,则认为是这个人,这个也是合理的。至此,也许没能讲清楚特征向量是什么,但我想对应于特征值大的特征向量有什么意义这一点也交代的差不多了。
3.2中交代了说了,这里想不到有什么需要补充的了。
理解了这几个问题之后,PCA简洁明了(这本来就是的)而又有深刻意义了。我突然发现原来看似简单的理论其实水深的很,o(︶︿︶)o !我看文章实在是不够认真,如果没有师兄提问,估计我也不会去想这个问题。再次感谢WF师兄。
 
[1]Eigenfaces for Recognition. Matthew Turk,Alex Pentland 1991
注:A的每一列是一张人脸,这里的人脸是每张人脸减去了所有人脸的平均后的人脸

Categories
不学无术 木有技术

第二章 啊哈!算法

原文见:

http://blog.csdn.net/silenough/article/details/7040028

 
A. 给定一个最多包含40亿个随机排列的32位整数的顺序文件,找出一个不在文件中的32位整数(在文件中至少缺少一个这样的数—为什么?)。在具有足够内存的情况下,如何解决该问题?如果有几个外部的“临时”文件可用,但是仅有几百字节的内村,又该如何解决该问题?
因为2^32 大于40亿,所以文件中至少缺失一个整数。我们从表示每个整数的32位的视角来考虑二分搜索,算法的第一趟(最多)读取40亿个输入整数,并把起始位为0的整数写入一个顺序文件,把起始位为1的整数写入另一个顺序文件。这两个文件中,有一个文件最多包含20亿个整数,接下来将该文件用作当前输入并重复探测过程,但这次探测的是第二个位。如果原始的输入文件包含n个元素,那么第一趟将读取n个整数,第二趟最多读取n/2个整数,第三趟最多读取n/4个整数,依次类推,所以总的运行时间正比于n。
如果内存足够,采用位图技术,通过排序文件并扫描,也能够找到缺失的整数,但是这样做会导致运行时间正比于nlog(n).
1. 考虑查找给定输入单词的所有变位词的问题。仅给定单词和字典的情况下,如何解决该问题?如果有一些时间和空间可以在响应任何查询之前预先处理字典,又会如何?
首先计算给定单词的标识,若果不允许预处理,那么久只能顺序读取整个文件,计算每个单词的标识,并于给定单词的标识进行比较。
 

  1. //压缩一个单词,形成其标识,设定单词中相同字母不会超过99个
  2. void compress(char * pWord, int len, char * pFlag)
  3. {
  4.     sort(pWord, pWord+len); //对单词进行排序
  5.     int i = 0;
  6.     int nCount;            //计数重复字母的个数
  7.     char chCount[3];       //存放整数到字符的转换值,整数最大为99
  8.     while (*pWord != ‘’)
  9.     {
  10.         char chTemp = *pWord;
  11.         char *pTemp = pWord + 1;
  12.         nCount = 1;
  13.         while (true)
  14.         {
  15.             if (chTemp == *pTemp++)
  16.             {
  17.                 ++nCount;
  18.             }
  19.             else
  20.             {
  21.                 *(pFlag + i++) = *pWord;
  22.             //  ++i;
  23.                 memset(chCount, ‘’, 3);
  24.                 _itoa(nCount, chCount, 10);
  25.                 if (nCount >= 10)
  26.                 {
  27.                     *(pFlag + i++) = *(chCount + 0);
  28.                 //  i++;
  29.                     *(pFlag + i++) = *(chCount + 1);
  30.                 //  ++i;
  31.                 }
  32.                 else
  33.                 {
  34.                     *(pFlag + i++) = *(chCount + 0);
  35.                 //  ++i;
  36.                 }
  37.                 pWord = pWord + nCount;
  38.                 break;
  39.             }
  40.         }
  41.     }
  42. }

如果允许进行预处理,我们可以在一个预先计算好的结构中执行二分查找,该结构中包含按标识排序的(标识,单词)对。
2. 给定包含4300 000 000 个32位整数的顺序文件,如何找出一个出现至少两次的整数?
 方法一:
二分搜索通过递归搜索包含半数以上整数的子区间来查找至少出现两次的单词。因为4300000000  > 2^32,所以必定存在重复的整数,搜索范围从[0, 2^32)开始,中间值mid为2^31,若区间[0, 2^31)内的整数个数大于2^31个,则调整搜索区间为[0, 2^31),反之则调整搜索区间为[2^31, 2^32),然后再对整个文件再遍历一遍,直到得到最后的结果。这样一共会有log2(n)次的搜索,每次遍历n个整数(每次都是完全遍历),总体的复杂度为o(nlog2(n))。
 

  1. #include <iostream>
  2. using namespace std;
  3. int pow2(int n)   //求2的n次幂
  4. {
  5.     int i;
  6.     int r = 1;
  7.     for (i = 0; i < n; i++)
  8.     {
  9.         r *=2;
  10.     }
  11.     return r;
  12. }
  13. int _tmain(int argc, _TCHAR* argv[])
  14. {
  15.     int arr[] = {4,2,5,1,6,3,8,0,7,6,11,12,14,9,15,10,13};
  16.     int len = sizeof(arr) / sizeof(int);
  17.     int nCount = 0;
  18.     int bit = 4;
  19.     int low = 0;
  20.     int high = pow2(bit);
  21.     int mid = (low +high) / 2;
  22.     int i;
  23.     while (low <= high )
  24.     {
  25.         mid = (low + high) / 2;  //取中间值
  26.         nCount = 0;
  27.         for (i = 0; i < len; i++)    //计数[low, mid)范围内整数的个数
  28.         {
  29.             if (arr[i] < mid && arr[i] >= low)
  30.             {
  31.                 ++nCount;
  32.             }
  33.         }  //end for
  34.         if (nCount == 0)     //若nCount为0,则mid即为重复的整数
  35.         {
  36.             cout << mid<<endl;
  37.             break;
  38.         }
  39.         else
  40.         {
  41.             if (nCount > (mid – low))  //若大于mid与low的差值,
  42.             {                          //表明重复的整数落在区间[low, mid)
  43.                 high = mid;        //缩小区间
  44.             }
  45.             else
  46.             {
  47.                 low = mid;
  48.             }  //end if
  49.         } //end if () else()
  50.     }  //end while
  51. }

 
   方法二:
 

  1. #include <iostream>
  2. #include <algorithm>
  3. using namespace std;
  4. int _tmain(int argc, _TCHAR* argv[])
  5. {
  6.     int arr[] = {4,2,5,1,7,3,8,0,7,6,11,12,14,9,15,10,13};
  7.     int len = sizeof(arr) / sizeof(int);
  8.     sort(arr, arr + len);   //先进行排序
  9.     int i;
  10.     int increase = arr[0];
  11.     for (i = 0; i < len; i++)
  12.     {
  13.         if (arr[i] > (i + increase))
  14.         {
  15.             increase += (arr[i] – i – increase);
  16.             continue;
  17.         }
  18.         if (arr[i] < (i + increase))
  19.         {
  20.             cout << arr[i] << endl;
  21.             break;
  22.         }
  23.     }
  24. }

3. 前面涉及了两个需要精巧代码来实现的向量旋转算法,将其分别作为独立的程序实现。在每个程序中,i和n的最大公约数如何实现?
采用辗转相除法求两个整数的最大公约数。
 

  1. int gcd(int a, int b)
  2. {
  3.     int temp;
  4.     if (a < b)  //使a始终为最大数
  5.     {
  6.         temp = a;
  7.         a = b;
  8.         b = temp;
  9.     }
  10.     while (b != 0)
  11.     {
  12.         temp = a % b;
  13.         a = b;
  14.         b = temp;
  15.     }
  16.     return a;
  17. }

 
 
方法一:海豚算法
 

  1. void Shifting(char * pArry, int rotdistance, int len)
  2. {
  3.     int i, j;
  4.     char temp;
  5.     int igcd = gcd(rotdistance, len);
  6.     for (i = 0; i < igcd; i++)
  7.     {
  8.         temp = pArry[i];
  9.         j = i;
  10.         for (; 😉
  11.         {
  12.             int k = j + rotdistance;
  13.             k %= len;
  14.             if ( k == i)
  15.             {
  16.                 break;
  17.             }
  18.             pArry[j] = pArry[k];
  19.             j = k;
  20.         }
  21.         pArry[j] = temp;
  22.     }
  23. }

方法二:块交换算法
 

  1. #include <iostream>
  2. using namespace std;
  3. //交换pArry[a…a+m-1]和pArry[b…b+m-1]
  4. void myswap(char *pArry, int a, int b, int m)
  5. {
  6.     char temp;
  7.     for (int i = 0; i < m; i++)
  8.     {
  9.         temp = pArry[a + i];
  10.         pArry[a + i] = pArry[b + i];
  11.         pArry[b + i] = temp;
  12.     }
  13. }
  14. void Shifting(char * pArry, int rotdistance, int len)
  15. {
  16.     if (rotdistance == 0 || rotdistance == len)
  17.     {
  18.         return;
  19.     }
  20.     int i, j, p;
  21.     i = p = rotdistance;
  22.     j = len – p;
  23.     while (i != j)
  24.     {
  25.         if (i > j)
  26.         {
  27.             myswap(pArry, p – i, p, j);
  28.             i -= j;
  29.         }
  30.         else
  31.         {
  32.             myswap(pArry, p – i, p + j – i, i);
  33.             j -= i;
  34.         }
  35.     }
  36.     myswap(pArry, p – i, p, i);
  37. }
  38. int _tmain(int argc, _TCHAR* argv[])
  39. {
  40.     char arry[] = “abcdefghijklmn”;
  41.     int len = strlen(arry);
  42.     Shifting(arry, 10, len);
  43.     return 0;
  44. }

方法三:求逆算法
根据矩阵的转置理论,对于矩阵AB,要得到BA,则分别求A和B的转置A’, B’,然后对(A’B’)转置,即(A’B’)’ = BA。同理,可以得到另一种一维向量向左旋转的算法。将要被旋转的向量x看做两部分ab,这里a代表x中的前rotdistance个元素。首先对a部分进行反转,再对b部分进行反转,最后对整个向量x进行反转即可。
对于字符串“abcdefgh”, rotdistance = 3, len = 8:
reverse(1, rotdistance);          //cbadefgh
reverse(rotdistance+1, len);  //cbahgfed
reverse(1, len);                       //defghabc
 

  1. #include <iostream>
  2. using namespace std;
  3. //对字符串中第i个字符到第j个字符进行反转,i、j>=1
  4. void MyReverse(char * pArry, int i, int j)
  5. {
  6.     int front = i;
  7.     int tail = j;
  8.     char temp;
  9.     while (front != tail && front < tail)
  10.     {
  11.         temp = pArry[front – 1];
  12.         pArry[front – 1] = pArry[tail – 1];
  13.         pArry[tail – 1] = temp;
  14.         ++front;
  15.         –tail;
  16.     }
  17. }
  18. //将字符串左旋转rotdistance个字符
  19. void Shifting(char * pArry, int rotdistance, int len)
  20. {
  21.     if (rotdistance == 0 || rotdistance == len)
  22.     {
  23.         return;
  24.     }
  25.     MyReverse(pArry, 1, rotdistance);
  26.     MyReverse(pArry, rotdistance + 1, len);
  27.     MyReverse(pArry, 1, len);
  28. }
  29. int _tmain(int argc, _TCHAR* argv[])
  30. {
  31.     char arry[] = “abcdefgh”;
  32.     int len = strlen(arry);
  33.     Shifting(arry, 5, len);
  34.     cout << arry << endl;
  35.     return 0;
  36. }

 
5. 向量旋转函数将向量ab变为ba。如何将向量abc变成cba?(这个交换非相邻内存块的问题进行了建模)
可以将bc看做一个整体,然后运用向量旋转算法,得到bca。然后对bc运用向量旋转算法,得到cb。最后变换后的向量为即cba.
 

  1. //交换pArry[a…a+m-1]和pArry[b…b+m-1]
  2. void myswap(char *pArry, int a, int b, int m)
  3. {
  4.     char temp;
  5.     for (int i = 0; i < m; i++)
  6.     {
  7.         temp = pArry[a + i];
  8.         pArry[a + i] = pArry[b + i];
  9.         pArry[b + i] = temp;
  10.     }
  11. }

 
 

  1. //对向量pArry中起始于ibegainPos位置的irotdistance-ibegainPos个元素
  2. //与起始于ibegainPos+irotdistance位置到位置iendPos之前的元素进行交换
  3. void SuccessiveSwap(char * pArry, int ibegainPos, int irotdistance, int iendPos)
  4. {
  5.     int i, j;
  6.     i = irotdistance – ibegainPos;
  7.     j = iendPos – irotdistance;
  8.     while (i != j)
  9.     {
  10.         if (i > j)
  11.         {
  12.             myswap(pArry, irotdistance – i, irotdistance, j);
  13.             i -= j;
  14.         }
  15.         else
  16.         {
  17.             myswap(pArry, irotdistance – i, irotdistance + j – i, i);
  18.             j -= i;
  19.         }
  20.     }
  21.     myswap(pArry, irotdistance – i, irotdistance, i);
  22. }

6. 如何实现一个以名字的按键编码为参数,并返回所有可能的匹配名字的函数
用按键编码标识每一个名字,并根据标识排序,然后顺序读取排序后的文件并输出具有不同名字的相同标识。为了检索出给定按钮编码的名字,可以使用一种包含编码标识和其他数据的结构。然后对该结构排序,使用二分搜索查询按键编码。
7. 转置一个存储在磁带上的4000×4000的矩阵(每条记录的格式相同,为数十个字节)。如何将运行的时间减少到半个小时?
给每条记录插入列号和行号,然后调用系统的磁带排序程序先按列排序再按行排序,最后使用另一个程序删除列号和行号。
8. 给定一个n元实数集合、一个实数t和一个整数k,如何快速确定是否存在一个k元子集,其元素之和不超过t?
对n元实数集合先进行排序,然后计算前k个元素的和既可以确定是否存在这样一个子集。若采用快速排序,时间复杂度为nlog10(n)。
9. 顺序搜素和二分搜索代表了搜索时间和预处理时间之间的折中。处理一个n元表格时,需要执行多少次二分搜索才能弥补对表进行排序所消耗的预处理时间?
对于顺序搜索,搜索k次的时间复杂度为O(kn);若采用二分搜索则需要先排序,则二分搜索的时间复杂度为O(nlog10(n)+log2(n))
变位词程序的实现
 

  1. //对从文件中读入的每个单词调用qsort库函数排序,输出到新的文件中
  2. void mysign(FILE * pFile1, FILE * pFile2)
  3. {
  4.     char word[20];
  5.     char sig[20];
  6.     pFile1 = fopen(“..file1.txt”, “r”);
  7.     if ( NULL == pFile1)
  8.     {
  9.         cout << “Open file1 error!” << endl;
  10.         return ;
  11.     }
  12.     pFile2 = fopen(“..file2.txt”, “w”);
  13.     if (NULL == pFile2)
  14.     {
  15.         cout << “Open file2 error!” << endl;
  16.         return ;
  17.     }
  18.     while (!feof(pFile1))
  19.     {
  20.         memset(word, ‘’, 20);
  21.         memset(sig, ‘’, 20);
  22.         fscanf(pFile1, “%s”, word);
  23.         strncpy(sig, word, strlen(word));
  24.         qsort(sig, strlen(sig), sizeof(char), charcomp);
  25.         fprintf(pFile2, “%s   %sn”, sig, word);
  26.     }
  27.     fclose(pFile1);
  28.     fclose(pFile2);
  29. }


 
 
使用sort程序将具有相同标识的单词归拢到一起,形成一个新的文件。最后调用squash()函数将具有相同变位词的单词在同一行输出。
 

  1. //将具有相同变位词的单词在同一行输出
  2. void squash()
  3. {
  4.     FILE * pFile3 = fopen(“..file3.txt”, “r”);
  5.     if ( NULL == pFile3)
  6.     {
  7.         cout << “Open file3 error!” << endl;
  8.         return ;
  9.     }
  10.     FILE * pFile4 = fopen(“..file4.txt”, “w”);
  11.     if (NULL == pFile4)
  12.     {
  13.         cout << “Open file4 error!” << endl;
  14.     }
  15.     char word[20];
  16.     char sig[20];
  17.     char oldsig[20];
  18.     int linenum = 0;
  19.     memset(oldsig, ‘’, 20);
  20.     while (!feof(pFile3))
  21.     {
  22.         memset(word, ‘’, 20);
  23.         memset(sig, ‘’, 20);
  24.         fscanf(pFile3, “%s %s”, sig, word);
  25.         if (strncmp(sig, oldsig, strlen(sig)) != 0 )
  26.         {
  27.             fprintf(pFile4, “n”);
  28.         }
  29.         strncpy(oldsig, sig, strlen(sig));
  30.         fprintf(pFile4, “%s “, word);
  31.     }
  32.     fclose(pFile3);
  33.     fclose(pFile4);
  34. }