learn.cpp

SMO工具箱

	}
	else if(kernelType ==1){
	    if(binaryFeature ==0){
			b1 =E1 +y1*(a1-lambda1)*power(1+dotProduct(example[e1],nonZeroFeature[e1],
				example[e1],nonZeroFeature[e1]),degree)+y2*(a2-lambda2)*power(1+dotProduct(example[e1],
				nonZeroFeature[e1],example[e2],nonZeroFeature[e2]),degree) +oldb;
			
			b2 =E2 +y1*(a1-lambda1)*power(1+dotProduct(example[e1],nonZeroFeature[e1],
				example[e2],nonZeroFeature[e2]),degree)+y2*(a2-lambda2)*power(1+dotProduct(example[e2],
				nonZeroFeature[e2],example[e2],nonZeroFeature[e2]),degree) +oldb;
		}
		else{
			b1 =E1 +y1*(a1 -lambda1)*power(1+nonZeroFeature[e1],degree)+y2*(a2-lambda2)*power(1+dotProduct(example[e1],nonZeroFeature[e1],
				example[e2],nonZeroFeature[e2]),degree) +oldb;
           	b2 =E2 +y1*(a1 -lambda1)*power(1+dotProduct(example[e1],nonZeroFeature[e1],
				example[e2],nonZeroFeature[e2]),degree)+y2*(a2-lambda2)*power(1+nonZeroFeature[e2],degree) +oldb;
		}
	}
	else if(kernelType ==2){
			if(binaryFeature ==0){
			interValue =dotProduct(example[e1],nonZeroFeature[e1],example[e1],
				nonZeroFeature[e1])-2*dotProduct(example[e1],nonZeroFeature[e1],
				example[e2],nonZeroFeature[e2])+dotProduct(example[e2],nonZeroFeature[e2],example[e2],
				nonZeroFeature[e2]);
		}
		else{
			interValue =nonZeroFeature[e1]-2*dotProduct(example[e1],nonZeroFeature[e1],
				example[e2],nonZeroFeature[e2])+nonZeroFeature[e2];
         }
		b1 =E1+y1*(a1-lambda1)+y2*(a2-lambda2)*exp(-interValue*rbfConstant)+oldb;
        b2 =E2+y1*(a1-lambda1)*exp(-interValue*rbfConstant)+y2*(a2-lambda2)+oldb;
	}
	if (fabs(b1-b2) <EPS)  //b==b1===b2
		b =b1;
	else if(!unBound1 &&!unBound2)
		b =(b1+b2)/2;
	else{
		if(unBound1)
			b=b1;
		if(unBound2)
			b=b2;
	}

	/********update the weight vector if kernel is linear******/

	if( kernelType==0){  //linear
		for(i =0;i<nonZeroFeature[e1];i++){
		if(binaryFeature ==0)
			 weight[example[e1][i]->id] +=y1*(a1-lambda1)*(example[e1][i]->value);
		else
             weight[example[e1][i]->id] +=y1*(a1-lambda1);
		}
	for(i =0;i<nonZeroFeature[e2];i++){
		if(binaryFeature ==0)
			weight[example[e2][i]->id] +=y2*(a2-lambda2)*(example[e2][i]->value);
		else
             weight[example[e2][i]->id] +=y2*(a2-lambda2);
	}
	}
	
	/****update the error cache and store a1 and a2.*******/
	for(i=0;i<numNonBound;i++){
		index = unBoundIndex[i];
		if(kernelType ==0)
			error[index] +=y1*(a1-lambda1)*dotProduct(example[e1],nonZeroFeature[e1],
			           example[index],nonZeroFeature[index]) + y2*(a2-lambda2)*dotProduct(example[e2],
					   nonZeroFeature[e2],example[index],nonZeroFeature[index])+oldb-b;
		else if(kernelType ==1)
            error[index] +=y1*(a1-lambda1)*power(1+dotProduct(example[e1],nonZeroFeature[e1],
			           example[index],nonZeroFeature[index]),degree) + y2*(a2-lambda2)*power(1+dotProduct(example[e2],
					   nonZeroFeature[e2],example[index],nonZeroFeature[index]),degree)+oldb-b;
		else if(kernelType ==2){
			if(binaryFeature ==0){
				interValue1 =dotProduct(example[e1],nonZeroFeature[e1],example[e1],nonZeroFeature[e1]) -2*dotProduct(
					 example[e1],nonZeroFeature[e1],example[index],nonZeroFeature[index])+dotProduct(example[index],nonZeroFeature[index],
					 example[index],nonZeroFeature[index]);

             interValue2 =dotProduct(example[e2],nonZeroFeature[e2],example[e2],nonZeroFeature[e2]) -2*dotProduct(
					 example[e2],nonZeroFeature[e2],example[index],nonZeroFeature[index])+dotProduct(example[index],nonZeroFeature[index],
					 example[index],nonZeroFeature[index]);
			}
			else{
				interValue1 = nonZeroFeature[e1] -2*dotProduct(example[e1],nonZeroFeature[e1],example[index],nonZeroFeature[index])
					          +nonZeroFeature[index];
				interValue2 = nonZeroFeature[e2] -2*dotProduct(example[e2],nonZeroFeature[e2],example[index],nonZeroFeature[index])
					          +nonZeroFeature[index];
			}
			error[index] +=y1*(a1-lambda1)*exp(-interValue1*rbfConstant)+y2*(a2-lambda2)*exp(-interValue2*rbfConstant) +oldb -b;
		}
	}
	lambda[e1] =a1;
	lambda[e2] =a2;
	
	/********calculate the error for e1 and e2*******/
	if(unBound1)
		error[e1]=0;
	if(unBound2)
		error[e2]=0;
	if(errorPtr>0)
		qsort2(errorCache,0,errorPtr,error);

	/*******update the nonBound set, the unBoundIndex set and the error cache index*****/
	if(!nonBound[e1] && unBound1){
		/*******************************
		    e1 was bound and is unbound after optimization.
			Add e1 to unboundIndex[],increment the number of nonbound and update nonBound[].
			Update error cache indexes and sort them in increasing order of error.
			**********************/
		unBoundPtr++;
		unBoundIndex[unBoundPtr] =e1;
		nonBound[e1] =1;
		numNonBound++;
		quicksort(unBoundIndex,0,unBoundPtr);
		errorPtr++;
		errorCache[errorPtr] =e1;
		if(errorPtr>0){
			qsort2(errorCache,0,errorPtr,error);
		}
	}
	else if(nonBound[e1] && !unBound1){
			/*******************************
		    e1 was nonbound and is bound after optimization.
			Remove e1 from unboundIndex[],decrement the number of nonbound and update nonBound[].
			Update error cache indexes and sort them in increasing order of error.
			**********************/
        findPos =binSearch(e1,unBoundIndex,numNonBound);
		   /*********************
		      Mark this previou e1 position in array unBoundIndex[] with 
			  a number greater than all possible indexes so that when we do sorting,this
			  number is at the end of range of the nonBound indexes.
		    ********************/
           unBoundIndex[findPos] = numExample+1;
		   quicksort(unBoundIndex,0,unBoundPtr);
		   unBoundPtr--;
		   numNonBound--;
		   nonBound[e1] =0;
		   if(errorCache[errorPtr] ==e1)
			   errorPtr--;
		   else{
			   error[e1] =error[errorCache[errorPtr]] +1;
			   qsort2(errorCache,0,errorPtr,error);
			   errorPtr--;
		   }
	}
    	if(!nonBound[e2] && unBound2){
		/*******************************
		    e2 was bound and is unbound after optimization.
			Add e2 to unboundIndex[],increment the number of nonbound and update nonBound[].
			Update error cache indexes and sort them in increasing order of error.
			**********************/
		unBoundPtr++;
		unBoundIndex[unBoundPtr] =e2;
		nonBound[e2] =1;
		numNonBound++;
		quicksort(unBoundIndex,0,unBoundPtr);
		errorPtr++;
		errorCache[errorPtr] =e2;
		if(errorPtr>0){
			qsort2(errorCache,0,errorPtr,error);
		}
	}
	else if(nonBound[e2] && !unBound2){
			/*******************************
		    e2 was nonbound and is bound after optimization.
			Remove e2 from unboundIndex[],decrement the number of nonbound and update nonBound[].
			Update error cache indexes and sort them in increasing order of error.
			**********************/
        findPos =binSearch(e2,unBoundIndex,numNonBound);
		   /*********************
		      Mark this previou e2 position in array unBoundIndex[] with 
			  a number greater than all possible indexes so that when we do sorting,this
			  number is at the end of range of the nonBound indexes.
		    ********************/
           unBoundIndex[findPos] = numExample+1;
		   quicksort(unBoundIndex,0,unBoundPtr);
		   unBoundPtr--;
		   numNonBound--;
		   nonBound[e2] =0;
		   //update error cache indices
		   if(errorCache[errorPtr] ==e2)
			   errorPtr--;
		   else{
			   error[e2] =error[errorCache[errorPtr]] +1;
			   qsort2(errorCache,0,errorPtr,error);
			   errorPtr--;
		   }
	}
	/**********iteration is the number of successful joint optimization********/

	iteration +=1;
	return 1;    
}

/*********examineExample()**********/
int examineExample(int e1)
{
	double r1;
	int found ,e2,i,index,y1;
	found =0;
	totalIteration++;
	y1 =target[e1];
	lambda1 =lambda[e1];
	if(nonBound[e1])
		E1 =error[e1];
	else 
		E1 =calculateError(e1);

	r1 =E1*y1;
	if( (r1<-TOL && lambda1<C) || (r1>TOL && lambda1>0)){
		//e1 violates KKT condition
		if(numNonBound >1){
			//this is the first hierachy of the second choice heuristic.
			if(E1>0){
				if(error[errorCache[0]]>=EPS)
					found =0;
				else{
					e2 =errorCache[0];
					found=1;
				}
			}
			else if(E1<0){
                if(error[errorCache[0]]<=EPS)
					found =0;
				else{
					e2 =errorCache[0];
					found=1;
				}
			}
			if(found)
				if(takeStep(e1,e2))
					return 1;
		}

		if(numNonBound>1){
			index =myrandom(numNonBound);
			e2 =unBoundIndex[index];
			for(i=0;i<numNonBound;i++){
				if(takeStep(e1,e2))
					return 1;
				index++;
				if(index ==numNonBound)
					index =0;
				e2=unBoundIndex[index];
			}
		}

		//this is third hierachy of the second choice heuristic.
		if( numNonZeroLambda>1){
			index = myrandom(numNonZeroLambda);
			e2 =nonZeroLambda[index];
			//only use bounded lambdas
			for(i=0;i<numNonZeroLambda;i++){
				if(nonBound[e1]==0){
					if(takeStep(e1,e2))
						return 1;
				}
				index++;
				if(index ==numNonZeroLambda)
					index =0;
				e2=nonZeroLambda[index];
			}
		}

		//this is the fourth hierachy of the second choice heuristic
		e2 = myrandom(numExample)+1;
		for(i =0;i<numExample;i++){
			//only use example with zero lambda
			if(lambda[i]<EPS){
				if(takeStep(e1,e2))
					return 1;
			}
			e2++;
			if (e2==numExample+1)
				e2=1;
		}
	}
	return 0;
} 

/**********public function***************
void startLearn(void)
Returns:Nothing.
Purpose:To find the support vectors from the training examples.
Notes:None.
**********/
void startLearn(void)
{
	int i;
	int numChanged=0;
	int examineAll=1;
	b =0; iteration =0; totalIteration =0;
    while(numChanged>0 ||examineAll){
		numChanged =0;
		if(examineAll){
			for(i =1;i<=numExample;i++)
				numChanged +=examineExample(i);
		}
		else {
			for(i =0;i<numNonBound;i++)
				numChanged +=examineExample(unBoundIndex[i]);
		}
		if(examineAll==1)
			examineAll=0;
		else if(numChanged ==0)
			examineAll =1;
	}
}