svm_struct_learn.c

SVMhmm: Learns a hidden Markov model from examples. Training examples (e.g. for part-of-speech taggi

      else {
	free_svector(diff);
      }

      if(struct_verbosity>=1)
	printf("(NumConst=%d, SV=%ld, CEps=%.4f, QPEps=%.4f)\n",cset.m,
	       svmModel->sv_num-1,ceps,svmModel->maxdiff);

      free_example(doc,0);
	
      rt_total+=MAX(get_runtime()-rt1,0);

  } while((ceps > sparm->epsilon) || 
	  finalize_iteration(ceps,cached_constraint,sample,sm,cset,alpha,sparm)
	 );
  

  if(struct_verbosity>=1) {
    /**** compute sum of slacks ****/
    /**** WARNING: If positivity constraints are used, then the
	  maximum slack id is larger than what is allocated
	  below ****/
    slacksum=0;
    if(sparm->slack_norm == 1) {
      for(j=0;j<cset.m;j++) 
	slacksum=MAX(slacksum,
		     cset.rhs[j]-classify_example(svmModel,cset.lhs[j]));
      }
    else if(sparm->slack_norm == 2) {
      exit(1);
    }
    alphasum=0;
    for(i=0; i<cset.m; i++)  
      alphasum+=alpha[i]*cset.rhs[i];
    modellength=model_length_s(svmModel,kparm);
    dualitygap=(0.5*modellength*modellength+sparm->C*(slacksum+ceps))
               -(alphasum-0.5*modellength*modellength);
    
    printf("Final epsilon on KKT-Conditions: %.5f\n",
	   MAX(svmModel->maxdiff,ceps));
    printf("Upper bound on duality gap: %.5f\n", dualitygap);
    printf("Dual objective value: dval=%.5f\n",
	    alphasum-0.5*modellength*modellength);
    printf("Total number of constraints in final working set: %i (of %i)\n",(int)cset.m,(int)totconstraints);
    printf("Number of iterations: %d\n",numIt);
    printf("Number of calls to 'find_most_violated_constraint': %ld\n",argmax_count);
    if(sparm->slack_norm == 1) {
      printf("Number of SV: %ld \n",svmModel->sv_num-1);
      printf("Norm of weight vector: |w|=%.5f\n",
	     model_length_s(svmModel,kparm));
    }
    else if(sparm->slack_norm == 2){ 
      printf("Number of SV: %ld (including %ld at upper bound)\n",
	     svmModel->sv_num-1,svmModel->at_upper_bound);
      printf("Norm of weight vector (including L2-loss): |w|=%.5f\n",
	     model_length_s(svmModel,kparm));
    }
    printf("Value of slack variable (on working set): xi=%.5f\n",slacksum);
    printf("Norm of longest difference vector: ||Psi(x,y)-Psi(x,ybar)||=%.5f\n",
	   length_of_longest_document_vector(cset.lhs,cset.m,kparm));
    printf("Runtime in cpu-seconds: %.2f (%.2f%% for QP, %.2f%% for kernel, %.2f%% for Argmax, %.2f%% for Psi, %.2f%% for init)\n",
	   rt_total/100.0, (100.0*rt_opt)/rt_total, (100.0*rt_kernel)/rt_total,
	   (100.0*rt_viol)/rt_total, (100.0*rt_psi)/rt_total, 
	   (100.0*rt_init)/rt_total);
  }
  if(ccache) {
    long cnum=0;
    CCACHEELEM *celem;
    for(i=0;i<n;i++) 
      for(celem=ccache->constlist[i];celem;celem=celem->next) 
	cnum++;
    printf("Final number of constraints in cache: %ld\n",cnum);
  }
  if(struct_verbosity>=4)
    printW(sm->w,sizePsi,n,lparm->svm_c);

  if(svmModel) {
    sm->svm_model=copy_model(svmModel);
    sm->w=sm->svm_model->lin_weights; /* short cut to weight vector */
  }

  print_struct_learning_stats(sample,sm,cset,alpha,sparm);

  if(ccache)    
    free_constraint_cache(ccache);
  for(i=0;i<n;i++)
    free_svector(fycache[i]);
  free(fycache);
  if(svmModel)
    free_model(svmModel,0);
  free(alpha); 
  free(alphahist); 
  free(cset.rhs); 
  for(i=0;i<cset.m;i++) 
    free_example(cset.lhs[i],1);
  free(cset.lhs);
  if(kparm->gram_matrix)
    free_matrix(kparm->gram_matrix);
}



void remove_inactive_constraints(CONSTSET *cset, double *alpha, 
			         long currentiter, long *alphahist, 
				 long mininactive)
     /* removes the constraints from cset (and alpha) for which
	alphahist indicates that they have not been active for at
	least mininactive iterations */

{  
  long i,m;
  
  m=0;
  for(i=0;i<cset->m;i++) {
    if((alphahist[i]<0) || ((currentiter-alphahist[i]) < mininactive)) {
      /* keep constraints that are marked as -1 or which have recently
         been active */
      cset->lhs[m]=cset->lhs[i];      
      cset->lhs[m]->docnum=m;
      cset->rhs[m]=cset->rhs[i];
      alpha[m]=alpha[i];
      alphahist[m]=alphahist[i];
      m++;
    }
    else {
      free_example(cset->lhs[i],1);
    }
  }
  if(cset->m != m) {
    cset->m=m;
    cset->lhs=(DOC **)realloc(cset->lhs,sizeof(DOC *)*cset->m);
    cset->rhs=(double *)realloc(cset->rhs,sizeof(double)*cset->m);
    /* alpha=realloc(alpha,sizeof(double)*cset->m); */
    /* alphahist=realloc(alphahist,sizeof(long)*cset->m); */
  }
}


MATRIX *init_kernel_matrix(CONSTSET *cset, KERNEL_PARM *kparm) 
     /* assigns a kernelid to each constraint in cset and creates the
	corresponding kernel matrix. */
{
  int i,j;
  CFLOAT kval;
  MATRIX *matrix;

  /* assign kernel id to each new constraint */
  for(i=0;i<cset->m;i++) 
    cset->lhs[i]->kernelid=i;

  /* allocate kernel matrix as necessary */
  matrix=create_matrix(i+50,i+50);

  for(j=0;j<cset->m;j++) {
    for(i=j;i<cset->m;i++) {
      kval=kernel(kparm,cset->lhs[j],cset->lhs[i]);
      matrix->element[j][i]=kval;
      matrix->element[i][j]=kval;
    }
  }
  return(matrix);
}

MATRIX *update_kernel_matrix(MATRIX *matrix, int newpos, CONSTSET *cset, 
			     KERNEL_PARM *kparm) 
     /* assigns new kernelid to constraint in position newpos and
	fills the corresponding part of the kernel matrix */
{
  int i,maxkernelid=0,newid;
  CFLOAT kval;
  double *used;

  /* find free kernelid to assign to new constraint */
  for(i=0;i<cset->m;i++) 
    if(i != newpos)
      maxkernelid=MAX(maxkernelid,cset->lhs[i]->kernelid);
  used=create_nvector(maxkernelid+2);
  clear_nvector(used,maxkernelid+2);
  for(i=0;i<cset->m;i++) 
    if(i != newpos)
      used[cset->lhs[i]->kernelid]=1;
  for(newid=0;used[newid];newid++);
  free_nvector(used);
  cset->lhs[newpos]->kernelid=newid;

  /* extend kernel matrix if necessary */
  maxkernelid=MAX(maxkernelid,newid);
  if((!matrix) || (maxkernelid>=matrix->m))
    matrix=realloc_matrix(matrix,maxkernelid+50,maxkernelid+50);

  for(i=0;i<cset->m;i++) {
    kval=kernel(kparm,cset->lhs[newpos],cset->lhs[i]);
    matrix->element[newid][cset->lhs[i]->kernelid]=kval;
    matrix->element[cset->lhs[i]->kernelid][newid]=kval;
  }
  return(matrix);
}

CCACHE *create_constraint_cache(SAMPLE sample, STRUCT_LEARN_PARM *sparm)
     /* create new constraint cache for training set */
{
  long        n=sample.n;
  EXAMPLE     *ex=sample.examples;
  CCACHE      *ccache;
  int         i;

  ccache=(CCACHE *)malloc(sizeof(CCACHE));
  ccache->n=n;
  ccache->constlist=(CCACHEELEM **)malloc(sizeof(CCACHEELEM *)*n);
  for(i=0;i<n;i++) { 
    /* add constraint for ybar=y to cache */
    ccache->constlist[i]=(CCACHEELEM *)malloc(sizeof(CCACHEELEM));
    ccache->constlist[i]->fydelta=create_svector_n(NULL,0,"",1);
    ccache->constlist[i]->rhs=loss(ex[i].y,ex[i].y,sparm)/n;
    ccache->constlist[i]->viol=0;
    ccache->constlist[i]->next=NULL;
  }
  return(ccache);
}

void free_constraint_cache(CCACHE *ccache)
     /* frees all memory allocated for constraint cache */
{
  CCACHEELEM *celem,*next;
  int i;
  for(i=0; i<ccache->n; i++) {
    celem=ccache->constlist[i];
    while(celem) {
      free_svector(celem->fydelta);
      next=celem->next;
      free(celem);
      celem=next;
    }
  }
  free(ccache->constlist);
  free(ccache);
}

void add_constraint_to_constraint_cache(CCACHE *ccache, MODEL *svmModel, int exnum, SVECTOR *fydelta, double rhs, int maxconst)
     /* add new constraint fydelta*w>rhs for example exnum to cache,
	if it is more violated than the currently most violated
	constraint in cache. if this grows the number of constraint
	for this example beyond maxconst, then the most unused
	constraint is deleted. the funciton assumes that
	update_constraint_cache_for_model has been run. */
{
  double  viol;
  double  dist_ydelta;
  DOC     *doc_fydelta;
  CCACHEELEM *celem;
  int     cnum;

  doc_fydelta=create_example(1,0,1,1,fydelta);
  dist_ydelta=classify_example(svmModel,doc_fydelta);
  free_example(doc_fydelta,0);  
  viol=rhs-dist_ydelta;

  if((viol-0.000000000001) > ccache->constlist[exnum]->viol) {
    celem=ccache->constlist[exnum];
    ccache->constlist[exnum]=(CCACHEELEM *)malloc(sizeof(CCACHEELEM));
    ccache->constlist[exnum]->next=celem;
    ccache->constlist[exnum]->fydelta=fydelta;
    ccache->constlist[exnum]->rhs=rhs;
    ccache->constlist[exnum]->viol=viol;

    /* remove last constraint in list, if list is longer than maxconst */
    cnum=2;
    for(celem=ccache->constlist[exnum];celem && celem->next && celem->next->next;celem=celem->next)
      cnum++;
    if(cnum>maxconst) {
      free_svector(celem->next->fydelta);
      free(celem->next);
      celem->next=NULL;
    }
  }
  else {
    free_svector(fydelta);
  }
}

void update_constraint_cache_for_model(CCACHE *ccache, MODEL *svmModel)
     /* update the violation scores according to svmModel and find the
	most violated constraints for each example */
{ 
  int     i;
  double  progress=0,progress_old=0;
  double  maxviol=0;
  double  dist_ydelta;
  DOC     *doc_fydelta;
  CCACHEELEM *celem,*prev,*maxviol_celem,*maxviol_prev;

  for(i=0; i<ccache->n; i++) { /*** example loop ***/
	  
    progress+=10.0/ccache->n;
    if((struct_verbosity==1) && (((int)progress_old) != ((int)progress)))
      {printf("+");fflush(stdout); progress_old=progress;}
    if(struct_verbosity>=2)
      {printf("+"); fflush(stdout);}
    
    maxviol=0;
    prev=NULL;
    maxviol_celem=NULL;
    maxviol_prev=NULL;
    for(celem=ccache->constlist[i];celem;celem=celem->next) {
      doc_fydelta=create_example(1,0,1,1,celem->fydelta);
      dist_ydelta=classify_example(svmModel,doc_fydelta);
      free_example(doc_fydelta,0);
      celem->viol=celem->rhs-dist_ydelta;
      if((celem->viol > maxviol) || (!maxviol_celem)) {
	maxviol=celem->viol;
	maxviol_celem=celem;
	maxviol_prev=prev;
      }
      prev=celem;
    }
    if(maxviol_prev) { /* move max violated constraint to the top of list */
      maxviol_prev->next=maxviol_celem->next;
      maxviol_celem->next=ccache->constlist[i];
      ccache->constlist[i]=maxviol_celem;
    }
  }
}

double find_most_violated_joint_constraint_in_cache(CCACHE *ccache, SVECTOR **lhs, double *margin)
     /* constructs most violated joint constraint from cache. assumes
	that update_constraint_cache_for_model has been run. NOTE:
	this function returns only a shallow copy of the Psi vectors
	in lhs. So, do not use a deep free, otherwise the case becomes
	invalid. */
{
  double sumviol=0;
  int i;
  SVECTOR *fydelta;

  (*lhs)=NULL;
  (*margin)=0;

  /**** add all maximally violated fydelta to joint constraint ****/
  for(i=0; i<ccache->n; i++) { 
    fydelta=copy_svector_shallow(ccache->constlist[i]->fydelta);
    append_svector_list(fydelta,(*lhs));           /* add fydelta to lhs */
    (*lhs)=fydelta;
    (*margin)+=ccache->constlist[i]->rhs;   /* add loss to rhs */
    sumviol+=ccache->constlist[i]->viol;
  }

  return(sumviol);
}


double find_most_violated_joint_constraint_in_cache_old(int n, int cache_size, SVECTOR ***fydelta_cache, double **loss_cache, MODEL *svmModel, SVECTOR **lhs, double *margin)
{
  int     i,j;
  double  progress,progress_old;
  double  maxviol=0,sumviol,viol,lossval;
  double  dist_ydelta;
  SVECTOR *fydelta;
  DOC     *doc_fydelta;

  (*lhs)=NULL;
  (*margin)=0;
  sumviol=0;

  progress=0;
  progress_old=progress;

  for(i=0; i<n; i++) { /*** example loop ***/
	  
    progress+=10.0/n;
    if((struct_verbosity==1) && (((int)progress_old) != ((int)progress)))
      {printf("+");fflush(stdout); progress_old=progress;}
    if(struct_verbosity>=2)
      {printf("+"); fflush(stdout);}
    
    fydelta=NULL;
    lossval=0;
    for(j=0;j<cache_size;j++) {
      doc_fydelta=create_example(1,0,1,1,fydelta_cache[j][i]);
      dist_ydelta=classify_example(svmModel,doc_fydelta);
      free_example(doc_fydelta,0);

      viol=loss_cache[j][i]-dist_ydelta;
      if((viol > maxviol) || (!fydelta)) {
	fydelta=fydelta_cache[j][i];
	lossval=loss_cache[j][i];
	maxviol=viol;
      }
    }

    /**** add current fydelta to joint constraint ****/
    fydelta=copy_svector(fydelta);
    append_svector_list(fydelta,(*lhs));     /* add fydelta to lhs */
    (*lhs)=fydelta;
    (*margin)+=lossval;                      /* add loss to rhs */
    sumviol+=maxviol;
  }

  return(sumviol);
}